OpenMP #pragma omp parallel for slow

Question

我試圖改進我的C源代碼以並行執行。我有一個四核CPU，所以我認爲4是很多線程（一個用於CPU）來運行我的程序，而不是像順序代碼那樣進行優化。

但它不工作。我的代碼沒有OpenMP需要11分鐘才能執行，而並行代碼需要11分鐘。我報告我的所有來源，但並行代碼僅在getBasin()函數中。

#include <stdio.h> 
#include <stdlib.h> 
#include <math.h> 
#include <time.h> 
#include <omp.h> 

#define BLD "\x1B[1m" 
#define RED "\x1B[31m" 
#define GRN "\x1B[32m" 
#define RST "\x1B[0m" 

struct point{ 
    double x; 
    double v; 
    double t; 
    double E0; 
    double E; 
    double dE; 
} typedef point; 

struct parametri { 
    double gamma; 
    double epsilon; 
    double dt; 
    double t_star; 
    double t_basin; 
    double t_max; 
    double x0; 
    double v0; 
    int choice; 
    int alg[1]; 
    double dx; 
} typedef parametri; 


// Prototipi delle funzioni 
void Setup(); 
void getParams(parametri *pars); 
void getError(point *xv, parametri *pars, int *i, int k, int *n_passi, double *xn1, double *vn1); 
void getBasin(point *xv, parametri *pars, int *h, int *n_passi, double *xn1, double *vn1); 

double f(double x, double v, parametri *pars); 

void algEulero(point *xv, parametri *pars, double *xn1, double *vn1); 
void algEuleroCromer(point *xv, parametri *pars, double *xn1, double *vn1); 
void algPuntoDiMezzo(point *xv, parametri *pars, double *xn1, double *vn1); 
void algVerlet(point *xv, parametri *pars, double *xn1, double *vn1); 
void algRungeKutta2(point *xv, parametri *pars, double *xn1, double *vn1); 

int main(void) { 

    // Inizializzo il display per l'utente. Maggiori informazioni vedere funzione Setup(); 
    Setup(); 

    // Dichiaro e recupero i parametri richiesti per il corretto funzionamento del programma; 
    parametri pars; 
    getParams(&pars); 

    // Dichiaro e ricavo i parametri essenziali, annesse variabili di supporto; 
    int i, n_passi = pars.t_max/pars.dt; 
    double dt0, xn1, vn1, t_star = 5.; 
    point xv; 

    // Imposto i parametri iniziali; 
    xv.x = pars.x0; 
    xv.v = pars.v0; 
    xv.E0 = 0.5*(xv.v)*(xv.v) - (xv.x)*(xv.x)/8. + (xv.x)*(xv.x)*(xv.x)*(xv.x)/4.; 
    xv.E = xv.E0; 
    xv.dE = 0; 
    xv.t = 0; 
    pars.t_star = 5; 
    pars.t_basin = 60; 
    dt0 = pars.dt; 

    // Formato dell'output; 
    printf ("t\tx(t)\tv(t)\tE(t)\tdE\n"); 

    if ((pars.choice == 1) || (pars.choice == 3) || (pars.choice == 4)) { // L'utente ha deciso di affrontare il primo/terzo/quarto esercizio; 

     // Avverto l'utente che il processo sta per iniziare. Può risultare inutile su tempi brevi, ma efficace su tempi molto lunghi; 
     fprintf(stderr, "\nAvvio integrazione numerica... ");  

     if (pars.alg[0] == 1) { // L'utente ha selezionato l'algoritmo di Eulero; 
      for (i=0; i<n_passi; i++){ 
       algEulero(&xv, &pars, &xn1, &vn1); 
      } 
     } else if (pars.alg[0] == 2) { // L'utente ha selezionato l'algoritmo di EuleroCromer; 
      for (i=0; i<n_passi; i++){ 
        algEuleroCromer(&xv, &pars, &xn1, &vn1); 
      } 
     } else if (pars.alg[0] == 3) { // L'utente ha selezionato l'algoritmo di PuntoDiMezzo; 
      for (i=0; i<n_passi; i++){ 
       algPuntoDiMezzo(&xv, &pars, &xn1, &vn1); 
      } 
     } else if (pars.alg[0] == 4) { // L'utente ha selezionato l'algoritmo di Verlet; 
      for (i=0; i<n_passi; i++){ 
       algVerlet(&xv, &pars, &xn1, &vn1); 
      } 
     } else if (pars.alg[0] == 5) { // L'utente ha selezionato l'algoritmo di RungeKutta di ordine 2; 
      for (i=0; i<n_passi; i++) { 
       algRungeKutta2(&xv, &pars, &xn1, &vn1); 
      } 
     } 

     // Algoritmo terminato; 
     fprintf(stderr, "[%s%sDONE%s]\n\n", BLD, GRN, RST); 

    } else if (pars.choice == 2) { // L'utente ha deciso di affrontare il secondo esercizio; 

     // Seleziono il secondo algoritmo da confrontare 
     do { 
      fprintf(stderr, "> Selezionare il secondo algoritmo:\n"); 
      fprintf(stderr, "\t>> [1] Eulero\n"); 
      fprintf(stderr, "\t>> [2] EuleroCromer - RungeKutta (1 ordine)\n"); 
      fprintf(stderr, "\t>> [3] PuntoDiMezzo\n"); 
      fprintf(stderr, "\t>> [4] Verlet\n"); 
      fprintf(stderr, "\t>> [5] RungeKutta (2 ordine)\n"); 
      fprintf(stderr, "\t>> "); 
      scanf("%d", &pars.alg[1]); 
     } while ((pars.alg[1] <= 0)); 

     // Avverto l'utente che il processo sta per iniziare. Può risultare inutile su tempi brevi, ma efficace su tempi molto lunghi; 
     fprintf(stderr, "\nAvvio calcolo errori... "); 

     // Eseguo lo studio degli errori d'integrazione mediante il primo e secondo algoritmo scelto, rispettivamente nei file error_1.dat, error_2.dat; 
     getError(&xv, &pars, &i, 0, &n_passi, &xn1, &vn1); 

     // Resetto le variabili e richiamo l'algoritmo per calcolare gli errori; 
     pars.dt = dt0; 
     n_passi = pars.t_max/pars.dt; 
     xv.x = pars.x0; 
     xv.v = pars.v0; 
     xv.E = xv.E0; 
     xv.dE = 0; 
     xv.t = 0; 

     getError(&xv, &pars, &i, 1, &n_passi, &xn1, &vn1); 

     // Processo terminato; 
     fprintf(stderr, "\n\nStato: [%s%sDONE%s]\n\n", BLD, GRN, RST); 

    } else if (pars.choice == 5) { // L'utente ha deciso di affrontare il quinto esercizio; 

     // Avverto l'utente che il processo sta per iniziare. Può risultare inutile su tempi brevi, ma efficace su tempi molto lunghi; 
     fprintf(stderr, "\nAvvio calcolo griglia... "); 

     getBasin(&xv, &pars, &i, &n_passi, &xn1, &vn1); 

     // Processo terminato; 
     fprintf(stderr, "[%s%sDONE%s]\n\n", BLD, GRN, RST); 

    } else { // L'utente non ha selezionato un esercizio valido; 
     fprintf(stderr, "\n[%s%sFAILED%s] Esercizio non disponibile.\n", BLD, RED, RST); 
     exit(EXIT_FAILURE); 
    } 

    return 0; 
} 


void Setup() { 

    fprintf(stderr, "\nAnalisi numerica di un'equazione differenziale\n"); 
    fprintf(stderr, "==============================================\n\n"); 

} 

void getParams(parametri *pars) { 

    do { 
     fprintf(stderr, "> Inserire un valore per gamma  : "); 
     scanf("%lf", &pars->gamma); 
    } while (pars->gamma < 0); 

    do { 
     fprintf(stderr, "> Inserire un valore per epsilon : "); 
     scanf("%lf", &pars->epsilon); 
    } while (pars->epsilon < 0); 

    do { 
     fprintf(stderr, "> Inserire un valore per dt   : "); 
     scanf("%lf", &pars->dt); 
    } while (pars->dt < 0); 

    do { 
     fprintf(stderr, "> Inserire un valore per tmax t.c. :\n"); 
     fprintf(stderr, " >> (tmax > 5) per I eserc.\n"); 
     fprintf(stderr, " >> (tmax > 60) per V eserc.  : "); 
     scanf("%lf", &pars->t_max); 
    } while (pars->t_max < 0); 

    do { 
     fprintf(stderr, "> Inserire un valore per x(0)  : "); 
     scanf("%lf", &pars->x0); 
    } while (pars->x0 < -1); 

    do { 
     fprintf(stderr, "> Inserire un valore per v(0)  : "); 
     scanf("%lf", &pars->v0); 
    } while (pars->v0 < -1); 

    do { 
     fprintf(stderr, "> Selezionare l'esercizio richiesto :\n"); 
     fprintf(stderr, "\t>> [1] Esercizio 1\n"); 
     fprintf(stderr, "\t>> [2] Esercizio 2\n"); 
     fprintf(stderr, "\t>> [3] Esercizio 3\n"); 
     fprintf(stderr, "\t>> [4] Esercizio 4\n"); 
     fprintf(stderr, "\t>> [5] Esercizio 5\n\n"); 

     fprintf(stderr, "\t>> "); 
     scanf("%d", &pars->choice); 
    } while ((pars->choice <= 0)); 

    do { 
     fprintf(stderr, "> Selezionare l'algoritmo voluto :\n"); 
     fprintf(stderr, "\t>> [1] Eulero\n"); 
     fprintf(stderr, "\t>> [2] EuleroCromer - RungeKutta (1 ordine)\n"); 
     fprintf(stderr, "\t>> [3] PuntoDiMezzo\n"); 
     fprintf(stderr, "\t>> [4] Verlet\n"); 
     fprintf(stderr, "\t>> [5] RungeKutta (2 ordine)\n\n"); 
     fprintf(stderr, "\t>> "); 
     scanf("%d", &pars->alg[0]); 
    } while ((pars->alg[0] <= 0)); 

} 

void getError(point *xv, parametri *pars, int *i, int k, int *n_passi, double *xn1, double *vn1) { 
    void (*algF)(point *, parametri *, double *, double *); 
    int j, n = 0; 
    FILE *fp; 

    // Questo if controlla se il codice sta eseguendo lo studio degli errori per il primo o secondo algoritmo (pars->alg[k]); 
    if (k == 0) { 
     fp = fopen("error_1.dat", "w+"); 
    } else { 
     fp = fopen("error_2.dat", "w+"); 
    } 

    // Assegno il puntatore corretto a seconda dell'algoritmo scelto; 
    if (pars->alg[k] == 1) { 
     algF = algEulero; 
    } else if (pars->alg[k] == 2) { 
     algF = algEuleroCromer; 
    } else if (pars->alg[k] == 3) { 
     algF = algPuntoDiMezzo; 
    } else if (pars->alg[k] == 4) { 
     algF = algVerlet; 
    } else if (pars->alg[k] == 5) { 
     algF = algRungeKutta2; 
    } else { 
     fprintf(stderr, "\n[%s%sFAILED%s] E' stato selezionato un algoritmo non valido.\n", BLD, RED, RST); 
     exit(EXIT_FAILURE); 
    } 

    // Informo l'utente che il processo sta iniziando; 
    fprintf(stderr, "\n>> Avvio %d algoritmo... ", k+1); 

    // Formattazione dell'output del file contenente gli errori; 
    fprintf(fp, "dt\tE(t*)\tE(0)\tdE/E(0)\t# passi\ti\tt\n"); 

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

     for ((*i)=0; (*i)<(*n_passi); (*i)++){ 
      (*algF)(xv, pars, xn1, vn1); 

      if (((pars->t_star - pars->dt/2. <= xv->t) && (xv->t >= pars->t_star + pars->dt/2.)) && (n == 0)) { 
       fprintf(fp, "%+.14lf\t%+.14lf\t%+.14lf\t%+.14lf\t%+.14d\t%+.14d\t%+.14lf\n", pars->dt, xv->E, xv->E0, (xv->E - xv->E0)/xv->E0, (*n_passi), (*i), xv->t); 
       n = 1; 
      } 
     } 

     // Resetto le variabili per rilanciare l'algoritmo con dt/2^j 
     n = 0; 
     xv->t = 0; 
     xv->x = pars->x0; 
     xv->v = pars->v0; 
     pars->dt = pars->dt/2.; 
     (*n_passi) = pars->t_max/pars->dt; 
     (*xn1) = 0; 
     (*vn1) = 0; 
     xv->E = xv->E0; 
    } 

    fclose(fp); 

    fprintf(stderr, "[%s%sDONE%s]", BLD, GRN, RST); 
} 

void getBasin(point *xv, parametri *pars, int *h, int *n_passi, double *xn1, double *vn1) { 

    // Dichiaro e setto i parametri che delimitano la griglia; 
    point xv1; 
    pars->dx = 0.01; 
    xv1.x = -1; 
    xv1.v = -1; 

    // Dichiaro le variabili necessarie per il bacino d'attrazione; 
    int i, j, i_max = 2./pars->dx, j_max = 2./pars->dx; 
    void (*algF)(point *, parametri *, double *, double *); 
    FILE *fp = fopen("basin.dat", "w+"); 

    // Assegno il puntatore corretto a seconda dell'algoritmo scelto; 
    if (pars->alg[0] == 1) { 
     algF = algEulero; 
    } else if (pars->alg[0] == 2) { 
     algF = algEuleroCromer; 
    } else if (pars->alg[0] == 3) { 
     algF = algPuntoDiMezzo; 
    } else if (pars->alg[0] == 4) { 
     algF = algVerlet; 
    } else if (pars->alg[0] == 5) { 
     algF = algRungeKutta2; 
    } else { 
     fprintf(stderr, "\n[%s%sFAILED%s] E' stato selezionato un algoritmo non valido.\n", BLD, RED, RST); 
     exit(EXIT_FAILURE); 
    } 

    // Formattazione output file basin.dat; 
    fprintf(fp, "x(0)\tx'(0)\tx(t*)\tv(t*)\n"); 

    omp_set_num_threads(4); 

    #pragma omp parallel for 
    // Eseguo il for della griglia sull'asse x'; 
    for (j=0; j<=j_max; j++) {  

     // Eseguo il for della griglia sull'asse x; 
     for (i=0; i<=i_max; i++) { 
      fprintf(fp, "%lf\t%lf\t", xv1.x, xv1.v); 

      xv->x = xv1.x; 
      xv->v = xv1.v; 

      // Eseguo l'integrazione numerica 
      for ((*h)=0; (*h)<(*n_passi); (*h)++) { 
       (*algF)(xv, pars, xn1, vn1); 

       // Entro in t = t*, stampo v(t*) ed esco dal ciclo; 
       if ((pars->t_basin - pars->dt/2. <= xv->t) && (xv->t >= pars->t_basin + pars->dt/2.)) { 
        fprintf(fp, "%lf\t%lf\n", xv->x, xv->v); 
        break; 
       } 
      } 

      xv1.x += pars->dx; 
      xv->t = 0; 
      (*xn1) = 0; 
      (*vn1) = 0; 
     } 

     // Resetto la x e incremento la x'; 
     xv1.x = -1; 
     xv1.v += pars->dx; 
    } 

} 

double f(double x, double v, parametri *pars) { 
    return 0.25*x - x*x*x + (pars->gamma - pars->epsilon*(x*x))*v; 
} 

void algEulero(point *xv, parametri *pars, double *xn1, double *vn1){ 
    printf("%.14g\t%.14g\t%.14g\t%.14g\t%.14g\n", xv->t, xv->x, xv->v, xv->E, xv->dE); 
    *xn1 = xv->x + (xv->v)*(pars->dt); 
    *vn1 = xv->v + f(xv->x, xv->v, pars)*(pars->dt); 

    xv->E = 0.5*(xv->v)*(xv->v) - (xv->x)*(xv->x)/8. + (xv->x)*(xv->x)*(xv->x)*(xv->x)/4.; 
    xv->dE = fabs(xv->E-xv->E0)/xv->E0; 

    xv->t += (pars->dt); 
    xv->x = *xn1; 
    xv->v = *vn1; 
} 

void algEuleroCromer(point *xv, parametri *pars, double *xn1, double *vn1){ 
    printf("%.14g\t%.14g\t%.14g\t%.14g\t%.14g\n", xv->t, xv->x, xv->v, xv->E, xv->dE); 
    xv->v = xv->v + f(xv->x, xv->v, pars)*(pars->dt); 
    xv->x = xv->x + (xv->v)*(pars->dt); 

    xv->E = 0.5*(xv->v)*(xv->v) - (xv->x)*(xv->x)/8. + (xv->x)*(xv->x)*(xv->x)*(xv->x)/4.; 
    xv->dE = fabs(xv->E-xv->E0)/xv->E0; 

    xv->t += (pars->dt); 
} 

void algPuntoDiMezzo(point *xv, parametri *pars, double *xn1, double *vn1) { 
    printf("%.14g\t%.14g\t%.14g\t%.14g\t%.14g\n", xv->t, xv->x, xv->v, xv->E, xv->dE); 
    *vn1 = xv->v + f(xv->x, xv->v, pars)*(pars->dt); 
    xv->x = xv->x + (0.5*(xv->v + (*vn1)))*(pars->dt); 

    xv->E = 0.5*(xv->v)*(xv->v) - (xv->x)*(xv->x)/8. + (xv->x)*(xv->x)*(xv->x)*(xv->x)/4.; 
    xv->dE = fabs(xv->E-xv->E0)/xv->E0; 

    xv->t += (pars->dt); 
    xv->v = *vn1; 
} 

void algVerlet(point *xv, parametri *pars, double *xn1, double *vn1) { 
    printf("%.14g\t%.14g\t%.14g\t%.14g\t%.14g\n", xv->t, xv->x, xv->v, xv->E, xv->dE); 
    *xn1 = xv->x + xv->v*pars->dt + 0.5*(f(xv->x, xv->v, pars))*pars->dt*pars->dt; 
    *vn1 = xv->v + 0.5*(f(xv->x, xv->v, pars) + f((*xn1), xv->v, pars))*pars->dt; 

    xv->x = *xn1; 
    xv->v = *vn1; 

    xv->E = 0.5*(xv->v)*(xv->v) - (xv->x)*(xv->x)/8. + (xv->x)*(xv->x)*(xv->x)*(xv->x)/4.; 
    xv->dE = fabs(xv->E-xv->E0)/xv->E0; 

    xv->t += (pars->dt); 
} 


void algRungeKutta2(point *xv, parametri *pars, double *xn1, double *vn1) { 
     printf("%.14g\t%.14g\t%.14g\t%.14g\t%.14g\n", xv->t, xv->x, xv->v, xv->E, xv->dE); 
     *xn1 = xv->x + xv->v*pars->dt + 0.5*f(xv->x, xv->v, pars)*pars->dt*pars->dt; 
     *vn1 = xv->v + f(xv->x + 0.5*xv->v*pars->dt, xv->v + 0.5*f(xv->x, xv->v, pars)*pars->dt, pars)*pars->dt; 

     xv->x = *xn1; 
     xv->v = *vn1; 

     xv->E = 0.5*(xv->v)*(xv->v) - (xv->x)*(xv->x)/8. + (xv->x)*(xv->x)*(xv->x)*(xv->x)/4.; 
     xv->dE = fabs(xv->E-xv->E0)/xv->E0; 

     xv->t += (pars->dt); 
} 

------------------編輯-----------------

尊敬的吉爾， 我向你解釋我的程序的功能。這個程序的目標是用不同算法（algEulero，algEuleroCromer，algPuntoDiMezzo，algVerlet，algRungeKutta2）在數值上求解一個微分方程。它工作正常。物理方程爲d^2x/dt^2 = f（x，v，gamma，epsilon）。這個f（）就是你可以在我的代碼中找到的f（）。我的簡單C程序的「大」問題是他非常慢：當我選擇5'練習（pars.choice == 5）時，程序將完成：

1）計算在getBasin（））中，一個面積爲[-1，1] x [-1,1]的pars.dx增量; 2）在x和y上的每個增量將啓動算法，該算法用for的兩個起始條件（x，y）數值地求解微分方程。當算法達到asynthotic t *（pars.t_basin）時，他將在basin.dat中寫出x（t *）和v（t *）的輸出。 3）（x，y）將改變並再次轉到點（1）。

現在，可以使用以下參數進行測試：0.83，4，0.01，62，1，1，5，5。

我測試您的代碼，但並不比我的順序代碼更快（在11分鐘內， 例如）。爲了改善它：

1）盆地輸出的順序是微不足道的，因爲我將繪製依賴於3'和4'列值的圖像（x，y）在Gnuplot上）。 2）爲什麼在函數getBasin（）之前創建線程，而不是僅爲兩個for（）？創建線程。這不應該重複，每個線程，getBasin（）？

對不起，我的英語不好，並行編程，我試圖改進它在線閱讀教程。

Answer 1

那麼，你的代碼的主要問題是並行版本是（非常）錯誤的。您可以在parallel區域之外定義所有變量，但不要聲明它們中的任何一個private（甚至不包括循環索引）。

此外，由於您將getBasin()函數的所有參數作爲指針傳遞，因此之後將它們聲明爲私有變得更加棘手。

但是，對代碼的快速分析表明，儘管這些參數是作爲指針傳遞的，但實際上並不在意它們在退出例程時的值。此外，它似乎沒有數據依賴於你試圖平行的循環（雖然我沒有做一個全面的檢查）。我發現的唯一明確的依賴關係是關於您打開的輸出文件，這將在並行更新中保持順序順序並不容易...

因此，爲了解決你的代碼的並行化，這是我做的：

1. 通過值而不是通過引用傳遞的參數getBasin()。這會生成一些額外的數據副本，但由於這些數據需要聲明爲private，因此無論如何都需要副本。
2. 在parallel區域內呼叫getBasin()。事實上，在我看來，這樣做比較簡單，而不是在函數本身內部處理數據私有化和IO問題。
3. 打開每個線程的輸出文件的一個私人副本，名爲「basinXX.dat」，其中「XX」是當前線程的ID，左側填充0，這些文件將包含全局輸出當前線程。希望這會適合你。否則，處理打印的順序可能會有點棘手。
4. 在函數內部使用單個孤兒omp for指令來並行循環。

因此，代碼應該（希望）可以擴展得更好。但是，由於您沒有指明用於計算的輸入參數，因此我無法對其進行測試，既無法驗證正確性，也無法驗證性能。因此，它可以只慘遭失敗...

無論如何，這裏就是我想出了：

#include <stdio.h> 
#include <stdlib.h> 
#include <math.h> 
#include <time.h> 
#include <omp.h> 

#define BLD "\x1B[1m" 
#define RED "\x1B[31m" 
#define GRN "\x1B[32m" 
#define RST "\x1B[0m" 

struct point{ 
    double x; 
    double v; 
    double t; 
    double E0; 
    double E; 
    double dE; 
} typedef point; 

struct parametri { 
    double gamma; 
    double epsilon; 
    double dt; 
    double t_star; 
    double t_basin; 
    double t_max; 
    double x0; 
    double v0; 
    int choice; 
    int alg[1]; 
    double dx; 
} typedef parametri; 


// Prototipi delle funzioni 
void Setup(); 
void getParams(parametri *pars); 
void getError(point *xv, parametri *pars, int *i, int k, int *n_passi, double *xn1, double *vn1); 
void getBasin(point xv, parametri pars, int h, int n_passi, double xn1, double vn1); 

double f(double x, double v, parametri *pars); 

void algEulero(point *xv, parametri *pars, double *xn1, double *vn1); 
void algEuleroCromer(point *xv, parametri *pars, double *xn1, double *vn1); 
void algPuntoDiMezzo(point *xv, parametri *pars, double *xn1, double *vn1); 
void algVerlet(point *xv, parametri *pars, double *xn1, double *vn1); 
void algRungeKutta2(point *xv, parametri *pars, double *xn1, double *vn1); 

int main(void) { 

    // Inizializzo il display per l'utente. Maggiori informazioni vedere funzione Setup(); 
    Setup(); 

    // Dichiaro e recupero i parametri richiesti per il corretto funzionamento del programma; 
    parametri pars; 
    getParams(&pars); 

    // Dichiaro e ricavo i parametri essenziali, annesse variabili di supporto; 
    int i, n_passi = pars.t_max/pars.dt; 
    double dt0, xn1, vn1, t_star = 5.; 
    point xv; 

    // Imposto i parametri iniziali; 
    xv.x = pars.x0; 
    xv.v = pars.v0; 
    xv.E0 = 0.5*(xv.v)*(xv.v) - (xv.x)*(xv.x)/8. + (xv.x)*(xv.x)*(xv.x)*(xv.x)/4.; 
    xv.E = xv.E0; 
    xv.dE = 0; 
    xv.t = 0; 
    pars.t_star = 5; 
    pars.t_basin = 60; 
    dt0 = pars.dt; 

    // Formato dell'output; 
    printf ("t\tx(t)\tv(t)\tE(t)\tdE\n"); 

    if ((pars.choice == 1) || (pars.choice == 3) || (pars.choice == 4)) { // L'utente ha deciso di affrontare il primo/terzo/quarto esercizio; 

     // Avverto l'utente che il processo sta per iniziare. Può risultare inutile su tempi brevi, ma efficace su tempi molto lunghi; 
     fprintf(stderr, "\nAvvio integrazione numerica... ");  

     if (pars.alg[0] == 1) { // L'utente ha selezionato l'algoritmo di Eulero; 
      for (i=0; i<n_passi; i++){ 
       algEulero(&xv, &pars, &xn1, &vn1); 
      } 
     } else if (pars.alg[0] == 2) { // L'utente ha selezionato l'algoritmo di EuleroCromer; 
      for (i=0; i<n_passi; i++){ 
        algEuleroCromer(&xv, &pars, &xn1, &vn1); 
      } 
     } else if (pars.alg[0] == 3) { // L'utente ha selezionato l'algoritmo di PuntoDiMezzo; 
      for (i=0; i<n_passi; i++){ 
       algPuntoDiMezzo(&xv, &pars, &xn1, &vn1); 
      } 
     } else if (pars.alg[0] == 4) { // L'utente ha selezionato l'algoritmo di Verlet; 
      for (i=0; i<n_passi; i++){ 
       algVerlet(&xv, &pars, &xn1, &vn1); 
      } 
     } else if (pars.alg[0] == 5) { // L'utente ha selezionato l'algoritmo di RungeKutta di ordine 2; 
      for (i=0; i<n_passi; i++) { 
       algRungeKutta2(&xv, &pars, &xn1, &vn1); 
      } 
     } 

     // Algoritmo terminato; 
     fprintf(stderr, "[%s%sDONE%s]\n\n", BLD, GRN, RST); 

    } else if (pars.choice == 2) { // L'utente ha deciso di affrontare il secondo esercizio; 

     // Seleziono il secondo algoritmo da confrontare 
     do { 
      fprintf(stderr, "> Selezionare il secondo algoritmo:\n"); 
      fprintf(stderr, "\t>> [1] Eulero\n"); 
      fprintf(stderr, "\t>> [2] EuleroCromer - RungeKutta (1 ordine)\n"); 
      fprintf(stderr, "\t>> [3] PuntoDiMezzo\n"); 
      fprintf(stderr, "\t>> [4] Verlet\n"); 
      fprintf(stderr, "\t>> [5] RungeKutta (2 ordine)\n"); 
      fprintf(stderr, "\t>> "); 
      scanf("%d", &pars.alg[1]); 
     } while ((pars.alg[1] <= 0)); 

     // Avverto l'utente che il processo sta per iniziare. Può risultare inutile su tempi brevi, ma efficace su tempi molto lunghi; 
     fprintf(stderr, "\nAvvio calcolo errori... "); 

     // Eseguo lo studio degli errori d'integrazione mediante il primo e secondo algoritmo scelto, rispettivamente nei file error_1.dat, error_2.dat; 
     getError(&xv, &pars, &i, 0, &n_passi, &xn1, &vn1); 

     // Resetto le variabili e richiamo l'algoritmo per calcolare gli errori; 
     pars.dt = dt0; 
     n_passi = pars.t_max/pars.dt; 
     xv.x = pars.x0; 
     xv.v = pars.v0; 
     xv.E = xv.E0; 
     xv.dE = 0; 
     xv.t = 0; 

     getError(&xv, &pars, &i, 1, &n_passi, &xn1, &vn1); 

     // Processo terminato; 
     fprintf(stderr, "\n\nStato: [%s%sDONE%s]\n\n", BLD, GRN, RST); 

    } else if (pars.choice == 5) { // L'utente ha deciso di affrontare il quinto esercizio; 

     // Avverto l'utente che il processo sta per iniziare. Può risultare inutile su tempi brevi, ma efficace su tempi molto lunghi; 
     fprintf(stderr, "\nAvvio calcolo griglia... "); 

     #pragma omp parallel num_threads(4) 
     getBasin(xv, pars, i, n_passi, xn1, vn1); 

     // Processo terminato; 
     fprintf(stderr, "[%s%sDONE%s]\n\n", BLD, GRN, RST); 

    } else { // L'utente non ha selezionato un esercizio valido; 
     fprintf(stderr, "\n[%s%sFAILED%s] Esercizio non disponibile.\n", BLD, RED, RST); 
     exit(EXIT_FAILURE); 
    } 

    return 0; 
} 


void Setup() { 

    fprintf(stderr, "\nAnalisi numerica di un'equazione differenziale\n"); 
    fprintf(stderr, "==============================================\n\n"); 

} 

void getParams(parametri *pars) { 

    do { 
     fprintf(stderr, "> Inserire un valore per gamma  : "); 
     scanf("%lf", &pars->gamma); 
    } while (pars->gamma < 0); 

    do { 
     fprintf(stderr, "> Inserire un valore per epsilon : "); 
     scanf("%lf", &pars->epsilon); 
    } while (pars->epsilon < 0); 

    do { 
     fprintf(stderr, "> Inserire un valore per dt   : "); 
     scanf("%lf", &pars->dt); 
    } while (pars->dt < 0); 

    do { 
     fprintf(stderr, "> Inserire un valore per tmax t.c. :\n"); 
     fprintf(stderr, " >> (tmax > 5) per I eserc.\n"); 
     fprintf(stderr, " >> (tmax > 60) per V eserc.  : "); 
     scanf("%lf", &pars->t_max); 
    } while (pars->t_max < 0); 

    do { 
     fprintf(stderr, "> Inserire un valore per x(0)  : "); 
     scanf("%lf", &pars->x0); 
    } while (pars->x0 < -1); 

    do { 
     fprintf(stderr, "> Inserire un valore per v(0)  : "); 
     scanf("%lf", &pars->v0); 
    } while (pars->v0 < -1); 

    do { 
     fprintf(stderr, "> Selezionare l'esercizio richiesto :\n"); 
     fprintf(stderr, "\t>> [1] Esercizio 1\n"); 
     fprintf(stderr, "\t>> [2] Esercizio 2\n"); 
     fprintf(stderr, "\t>> [3] Esercizio 3\n"); 
     fprintf(stderr, "\t>> [4] Esercizio 4\n"); 
     fprintf(stderr, "\t>> [5] Esercizio 5\n\n"); 

     fprintf(stderr, "\t>> "); 
     scanf("%d", &pars->choice); 
    } while ((pars->choice <= 0)); 

    do { 
     fprintf(stderr, "> Selezionare l'algoritmo voluto :\n"); 
     fprintf(stderr, "\t>> [1] Eulero\n"); 
     fprintf(stderr, "\t>> [2] EuleroCromer - RungeKutta (1 ordine)\n"); 
     fprintf(stderr, "\t>> [3] PuntoDiMezzo\n"); 
     fprintf(stderr, "\t>> [4] Verlet\n"); 
     fprintf(stderr, "\t>> [5] RungeKutta (2 ordine)\n\n"); 
     fprintf(stderr, "\t>> "); 
     scanf("%d", &pars->alg[0]); 
    } while ((pars->alg[0] <= 0)); 

} 

void getError(point *xv, parametri *pars, int *i, int k, int *n_passi, double *xn1, double *vn1) { 
    void (*algF)(point *, parametri *, double *, double *); 
    int j, n = 0; 
    FILE *fp; 

    // Questo if controlla se il codice sta eseguendo lo studio degli errori per il primo o secondo algoritmo (pars->alg[k]); 
    if (k == 0) { 
     fp = fopen("error_1.dat", "w+"); 
    } else { 
     fp = fopen("error_2.dat", "w+"); 
    } 

    // Assegno il puntatore corretto a seconda dell'algoritmo scelto; 
    if (pars->alg[k] == 1) { 
     algF = algEulero; 
    } else if (pars->alg[k] == 2) { 
     algF = algEuleroCromer; 
    } else if (pars->alg[k] == 3) { 
     algF = algPuntoDiMezzo; 
    } else if (pars->alg[k] == 4) { 
     algF = algVerlet; 
    } else if (pars->alg[k] == 5) { 
     algF = algRungeKutta2; 
    } else { 
     fprintf(stderr, "\n[%s%sFAILED%s] E' stato selezionato un algoritmo non valido.\n", BLD, RED, RST); 
     exit(EXIT_FAILURE); 
    } 

    // Informo l'utente che il processo sta iniziando; 
    fprintf(stderr, "\n>> Avvio %d algoritmo... ", k+1); 

    // Formattazione dell'output del file contenente gli errori; 
    fprintf(fp, "dt\tE(t*)\tE(0)\tdE/E(0)\t# passi\ti\tt\n"); 

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

     for ((*i)=0; (*i)<(*n_passi); (*i)++){ 
      (*algF)(xv, pars, xn1, vn1); 

      if (((pars->t_star - pars->dt/2. <= xv->t) && (xv->t >= pars->t_star + pars->dt/2.)) && (n == 0)) { 
       fprintf(fp, "%+.14lf\t%+.14lf\t%+.14lf\t%+.14lf\t%+.14d\t%+.14d\t%+.14lf\n", pars->dt, xv->E, xv->E0, (xv->E - xv->E0)/xv->E0, (*n_passi), (*i), xv->t); 
       n = 1; 
      } 
     } 

     // Resetto le variabili per rilanciare l'algoritmo con dt/2^j 
     n = 0; 
     xv->t = 0; 
     xv->x = pars->x0; 
     xv->v = pars->v0; 
     pars->dt = pars->dt/2.; 
     (*n_passi) = pars->t_max/pars->dt; 
     (*xn1) = 0; 
     (*vn1) = 0; 
     xv->E = xv->E0; 
    } 

    fclose(fp); 

    fprintf(stderr, "[%s%sDONE%s]", BLD, GRN, RST); 
} 

void getBasin(point xv, parametri pars, int h, int n_passi, double xn1, double vn1) { 

    // Dichiaro e setto i parametri che delimitano la griglia; 
    point xv1; 
    pars.dx = 0.01; 
    xv1.x = -1; 
    xv1.v = -1; 

    // Dichiaro le variabili necessarie per il bacino d'attrazione; 
    int i, j, i_max = 2./pars.dx, j_max = 2./pars.dx; 
    void (*algF)(point *, parametri *, double *, double *); 
    char fname[13]; 
    sprintf(fname, "bassin%02d.dat", omp_get_thread_num()); 

    FILE *fp = fopen(fname, "w+"); 

    // Assegno il puntatore corretto a seconda dell'algoritmo scelto; 
    if (pars.alg[0] == 1) { 
     algF = algEulero; 
    } else if (pars.alg[0] == 2) { 
     algF = algEuleroCromer; 
    } else if (pars.alg[0] == 3) { 
     algF = algPuntoDiMezzo; 
    } else if (pars.alg[0] == 4) { 
     algF = algVerlet; 
    } else if (pars.alg[0] == 5) { 
     algF = algRungeKutta2; 
    } else { 
     fprintf(stderr, "\n[%s%sFAILED%s] E' stato selezionato un algoritmo non valido.\n", BLD, RED, RST); 
     exit(EXIT_FAILURE); 
    } 

    // Formattazione output file basin.dat; 
    fprintf(fp, "x(0)\tx'(0)\tx(t*)\tv(t*)\n"); 

    #pragma omp for 
    // Eseguo il for della griglia sull'asse x'; 
    for (j=0; j<=j_max; j++) {  

     // Eseguo il for della griglia sull'asse x; 
     for (i=0; i<=i_max; i++) { 
      fprintf(fp, "%lf\t%lf\t", xv1.x, xv1.v); 

      xv.x = xv1.x; 
      xv.v = xv1.v; 

      // Eseguo l'integrazione numerica 
      for (h=0; h<n_passi; h++) { 
       (*algF)(&xv, &pars, &xn1, &vn1); 

       // Entro in t = t*, stampo v(t*) ed esco dal ciclo; 
       if ((pars.t_basin - pars.dt/2. <= xv.t) && (xv.t >= pars.t_basin + pars.dt/2.)) { 
        fprintf(fp, "%lf\t%lf\n", xv.x, xv.v); 
        break; 
       } 
      } 

      xv1.x += pars.dx; 
      xv.t = 0; 
      xn1 = 0; 
      vn1 = 0; 
     } 

     // Resetto la x e incremento la x'; 
     xv1.x = -1; 
     xv1.v += pars.dx; 
    } 

} 

double f(double x, double v, parametri *pars) { 
    return 0.25*x - x*x*x + (pars->gamma - pars->epsilon*(x*x))*v; 
} 

void algEulero(point *xv, parametri *pars, double *xn1, double *vn1){ 
    //printf("%.14g\t%.14g\t%.14g\t%.14g\t%.14g\n", xv->t, xv->x, xv->v, xv->E, xv->dE); 
    *xn1 = xv->x + (xv->v)*(pars->dt); 
    *vn1 = xv->v + f(xv->x, xv->v, pars)*(pars->dt); 

    xv->E = 0.5*(xv->v)*(xv->v) - (xv->x)*(xv->x)/8. + (xv->x)*(xv->x)*(xv->x)*(xv->x)/4.; 
    xv->dE = fabs(xv->E-xv->E0)/xv->E0; 

    xv->t += (pars->dt); 
    xv->x = *xn1; 
    xv->v = *vn1; 
} 

void algEuleroCromer(point *xv, parametri *pars, double *xn1, double *vn1){ 
    //printf("%.14g\t%.14g\t%.14g\t%.14g\t%.14g\n", xv->t, xv->x, xv->v, xv->E, xv->dE); 
    xv->v = xv->v + f(xv->x, xv->v, pars)*(pars->dt); 
    xv->x = xv->x + (xv->v)*(pars->dt); 

    xv->E = 0.5*(xv->v)*(xv->v) - (xv->x)*(xv->x)/8. + (xv->x)*(xv->x)*(xv->x)*(xv->x)/4.; 
    xv->dE = fabs(xv->E-xv->E0)/xv->E0; 

    xv->t += (pars->dt); 
} 

void algPuntoDiMezzo(point *xv, parametri *pars, double *xn1, double *vn1) { 
    //printf("%.14g\t%.14g\t%.14g\t%.14g\t%.14g\n", xv->t, xv->x, xv->v, xv->E, xv->dE); 
    *vn1 = xv->v + f(xv->x, xv->v, pars)*(pars->dt); 
    xv->x = xv->x + (0.5*(xv->v + (*vn1)))*(pars->dt); 

    xv->E = 0.5*(xv->v)*(xv->v) - (xv->x)*(xv->x)/8. + (xv->x)*(xv->x)*(xv->x)*(xv->x)/4.; 
    xv->dE = fabs(xv->E-xv->E0)/xv->E0; 

    xv->t += (pars->dt); 
    xv->v = *vn1; 
} 

void algVerlet(point *xv, parametri *pars, double *xn1, double *vn1) { 
    //printf("%.14g\t%.14g\t%.14g\t%.14g\t%.14g\n", xv->t, xv->x, xv->v, xv->E, xv->dE); 
    *xn1 = xv->x + xv->v*pars->dt + 0.5*(f(xv->x, xv->v, pars))*pars->dt*pars->dt; 
    *vn1 = xv->v + 0.5*(f(xv->x, xv->v, pars) + f((*xn1), xv->v, pars))*pars->dt; 

    xv->x = *xn1; 
    xv->v = *vn1; 

    xv->E = 0.5*(xv->v)*(xv->v) - (xv->x)*(xv->x)/8. + (xv->x)*(xv->x)*(xv->x)*(xv->x)/4.; 
    xv->dE = fabs(xv->E-xv->E0)/xv->E0; 

    xv->t += (pars->dt); 
} 


void algRungeKutta2(point *xv, parametri *pars, double *xn1, double *vn1) { 
     //printf("%.14g\t%.14g\t%.14g\t%.14g\t%.14g\n", xv->t, xv->x, xv->v, xv->E, xv->dE); 
     *xn1 = xv->x + xv->v*pars->dt + 0.5*f(xv->x, xv->v, pars)*pars->dt*pars->dt; 
     *vn1 = xv->v + f(xv->x + 0.5*xv->v*pars->dt, xv->v + 0.5*f(xv->x, xv->v, pars)*pars->dt, pars)*pars->dt; 

     xv->x = *xn1; 
     xv->v = *vn1; 

     xv->E = 0.5*(xv->v)*(xv->v) - (xv->x)*(xv->x)/8. + (xv->x)*(xv->x)*(xv->x)*(xv->x)/4.; 
     xv->dE = fabs(xv->E-xv->E0)/xv->E0; 

     xv->t += (pars->dt); 
} 

---

編輯

有了新的輸入參數你給的，我是能夠測試的代碼，事實證明：

1. 我在我給你的代碼中有一個小錯誤（用於打開輸出文件，使用"fname"而不是fname）
2. 您的代碼花費大部分（全部）時間在printf()。

固定這兩個和編譯這樣後：

gcc -O3 -fopenmp -mtune=native -march=native -w bassin.c 

我2.12s

跑在我的雙核筆記本電腦的代碼更新。請自行嘗試。