try to fix openmp

1 files changed, 44 insertions, 51 deletions

diff --git a/src/lossdistrib.c b/src/lossdistrib.c
index 412437d..494fed7 100644
--- a/src/lossdistrib.c
+++ b/src/lossdistrib.c
@@ -491,6 +491,11 @@ void lossdistrib_joint_Z(const double *dp, const double *pp, const int *ndp,
     free(qmat);
 }
 
+static double* gshocked;
+static double* gshockedbar;
+static double* Sshocked;
+static double* Rshocked;
+
 void BCloss_recov_dist(const double *defaultprob, const int *dim1, const int *dim2,
                        const double *issuerweights, const double *recov,
                        const double *Z, const double *w, const int *n,
@@ -521,31 +526,26 @@ void BCloss_recov_dist(const double *defaultprob, const int *dim1, const int *di
     for(int t = 0; t < (*dim2); t++) {
         memset(Lw, 0, (*N) * (*n) * sizeof(double));
         memset(Rw, 0, (*N) * (*n) * sizeof(double));
-        static double* gshocked;
-        static double* Sshocked;
-        static double* Rshocked;
+
         #pragma omp threadprivate(gshocked, Sshocked, Rshocked)
         #pragma omp parallel
         {
             gshocked = malloc(3 * (*dim1) * sizeof(double));
             Sshocked = gshocked + *dim1;
             Rshocked = Sshocked + *dim1;
-        }
-        #pragma omp parallel for
-        for(int i = 0; i < *n; i++) {
-            for(int j = 0; j < (*dim1); j++) {
-                double g = defaultprob[j + (*dim1) * t];
-                gshocked[j] = shockprob(g, rho[j], Z[i], 0);
-                Sshocked[j] = shockseverity(1-recov[j], Z[i], rho[j], g);
-                Rshocked[j] = 1 - Sshocked[j];
+            #pragma omp for
+            for(int i = 0; i < *n; i++) {
+                for(int j = 0; j < (*dim1); j++) {
+                    double g = defaultprob[j + (*dim1) * t];
+                    gshocked[j] = shockprob(g, rho[j], Z[i], 0);
+                    Sshocked[j] = shockseverity(1-recov[j], Z[i], rho[j], g);
+                    Rshocked[j] = 1 - Sshocked[j];
+                }
+                lossdistrib(gshocked, dim1, issuerweights, Sshocked, N, N, defaultflag,
+                            Lw + i * (*N));
+                lossdistrib(gshocked, dim1, issuerweights, Rshocked, N, N, defaultflag,
+                            Rw + i * (*N));
             }
-            lossdistrib(gshocked, dim1, issuerweights, Sshocked, N, N, defaultflag,
-                        Lw + i * (*N));
-            lossdistrib(gshocked, dim1, issuerweights, Rshocked, N, N, defaultflag,
-                        Rw + i * (*N));
-        }
-        #pragma omp parallel
-        {
             free(gshocked);
         }
         cblas_dgemv(CblasColMajor, CblasNoTrans, *N, *n, alpha, Lw, *N, w, 1, beta, L + t * (*N), 1);
@@ -590,54 +590,47 @@ void BCloss_recov_trunc(const double *defaultprob, const int *dim1, const int *d
     double* Lw = malloc(T * (*n) * sizeof(double));
     const double alpha = 1;
     const double beta = 0;
-
     for(int t = 0; t < (*dim2); t++) {
         memset(Lw, 0, T * (*n) * sizeof(double));
-        static double* gshocked;
-        static double* gshockedbar;
-        static double* Sshocked;
-        static double* Rshocked;
 
-        #pragma omp threadprivate(gshocked, gshockedbar,Sshocked, Rshocked)
+        #pragma omp threadprivate(gshocked, gshockedbar, Sshocked, Rshocked)
         #pragma omp parallel
         {
             gshocked = malloc(4 * (*dim1) * sizeof(double));
             gshockedbar = gshocked + *dim1;
             Rshocked = gshockedbar + *dim1;
             Sshocked = Rshocked + *dim1;
-        }
-        #pragma omp for
-        for(int i = 0; i < *n; i++) {
 
-            double g;
-            for(int j = 0; j < (*dim1); j++){
-                g = defaultprob[j + (*dim1) * t];
-                gshocked[j] = shockprob(g, rho[j], Z[i], 0);
-                Sshocked[j] = shockseverity(1-recov[j], Z[i], rho[j], g);
-                gshockedbar[j] = 1 - gshocked[j];
-                Rshocked[j] = 1 - Sshocked[j];
-            }
+            #pragma omp for
+            for(int i = 0; i < *n; i++) {
+
+                double g;
+                for(int j = 0; j < (*dim1); j++){
+                    g = defaultprob[j + (*dim1) * t];
+                    gshocked[j] = shockprob(g, rho[j], Z[i], 0);
+                    Sshocked[j] = shockseverity(1-recov[j], Z[i], rho[j], g);
+                    gshockedbar[j] = 1 - gshocked[j];
+                    Rshocked[j] = 1 - Sshocked[j];
+                }
 
-            lossdistrib(gshocked, dim1, issuerweights, Sshocked,
-                        N, &T, defaultflag, Lw + i * T);
-            ER[i] = 0;
-            double Ktilde = *K - cblas_ddot(*dim1, issuerweights, 1, Sshocked, 1);
+                lossdistrib(gshocked, dim1, issuerweights, Sshocked,
+                            N, &T, defaultflag, Lw + i * T);
+                ER[i] = 0;
+                double Ktilde = *K - cblas_ddot(*dim1, issuerweights, 1, Sshocked, 1);
 
-            if(Ktilde > 0) {
-                int Ttilde = (int) floor(Ktilde * (*N))+1;
+                if(Ktilde > 0) {
+                    int Ttilde = (int) floor(Ktilde * (*N))+1;
 
-                double* Rw = calloc(2 * Ttilde, sizeof(double));
-                double* valR = Rw + Ttilde;
-                lossdistrib(gshockedbar, dim1, issuerweights, Rshocked,
-                            N, &Ttilde, defaultflag, Rw);
+                    double* Rw = calloc(2 * Ttilde, sizeof(double));
+                    double* valR = Rw + Ttilde;
+                    lossdistrib(gshockedbar, dim1, issuerweights, Rshocked,
+                                N, &Ttilde, defaultflag, Rw);
 
-                posK(Ttilde, Ktilde, lu, valR);
-                ER[i] = cblas_ddot(Ttilde, Rw, 1, valR, 1);
-                free(Rw);
+                    posK(Ttilde, Ktilde, lu, valR);
+                    ER[i] = cblas_ddot(Ttilde, Rw, 1, valR, 1);
+                    free(Rw);
+                }
             }
-        }
-        #pragma omp parallel
-        {
             free(gshocked);
         }
         cblas_dgemv(CblasColMajor, CblasTrans, T, *n, alpha, Lw, T, valL, 1, beta, EL, 1);