[Getfem-commits] (no subject)

[Andriy Andreykiv]

branch: upgrading_getfem_to_cpp14
commit 701030ad2b8417e02172cc468b36d4ae50fcbebd
Author: aa <address@hidden>
Date:   Mon Feb 18 15:37:46 2019 +0100

    replacing boost::thread_specific_ptr with thread_local via a macro
---
 src/bgeot_geometric_trans.cc    |  8 ++--
 src/getfem/bgeot_tensor.h       |  8 ++--
 src/getfem/dal_basic.h          |  4 +-
 src/getfem/getfem_omp.h         | 30 +------------
 src/getfem_fem.cc               | 98 ++++++++++++++++++++---------------------
 src/getfem_integration.cc       | 24 +++++-----
 src/getfem_linearized_plates.cc |  4 +-
 src/getfem_mat_elem.cc          |  8 ++--
 8 files changed, 79 insertions(+), 105 deletions(-)

diff --git a/src/bgeot_geometric_trans.cc b/src/bgeot_geometric_trans.cc
index 0b1cc6a..7e89d4f 100644
--- a/src/bgeot_geometric_trans.cc
+++ b/src/bgeot_geometric_trans.cc
@@ -29,22 +29,22 @@
 namespace bgeot {
 
   std::vector<scalar_type>& __aux1(){
-    DEFINE_STATIC_THREAD_LOCAL(std::vector<scalar_type>, v);
+    THREAD_SAFE_STATIC std::vector<scalar_type> v;
     return v;
   }
 
   std::vector<scalar_type>& __aux2(){
-    DEFINE_STATIC_THREAD_LOCAL(std::vector<scalar_type>, v);
+    THREAD_SAFE_STATIC std::vector<scalar_type> v;
     return v;
   }
 
   std::vector<scalar_type>& __aux3(){
-    DEFINE_STATIC_THREAD_LOCAL(std::vector<scalar_type>, v);
+    THREAD_SAFE_STATIC std::vector<scalar_type> v;
     return v;
   }
 
   std::vector<long>& __ipvt_aux(){
-    DEFINE_STATIC_THREAD_LOCAL(std::vector<long>, vi);
+    THREAD_SAFE_STATIC std::vector<long> vi;
     return vi;
   }
 
diff --git a/src/getfem/bgeot_tensor.h b/src/getfem/bgeot_tensor.h
index 43d779d..b155d2b 100644
--- a/src/getfem/bgeot_tensor.h
+++ b/src/getfem/bgeot_tensor.h
@@ -333,8 +333,8 @@ namespace bgeot {
     /* reduction du tenseur t par son indice ni et la matrice          */
     /* transposee de m.                                                */
 
-    DEFINE_STATIC_THREAD_LOCAL(std::vector<T>, tmp);
-    DEFINE_STATIC_THREAD_LOCAL(multi_index, mi);
+    THREAD_SAFE_STATIC std::vector<T> tmp;
+    THREAD_SAFE_STATIC multi_index mi;
 
     mi = t.sizes();
     size_type dimt = mi[ni], dim = m.nrows();
@@ -400,8 +400,8 @@ namespace bgeot {
   template<class T> void tensor<T>::mat_reduction
   (const tensor &t, const gmm::dense_matrix<T> &m, int ni) {
     /* reduction du tenseur t par son indice ni et la matrice m.       */
-    DEFINE_STATIC_THREAD_LOCAL(std::vector<T>, tmp);
-    DEFINE_STATIC_THREAD_LOCAL(multi_index, mi);
+    THREAD_SAFE_STATIC std::vector<T> tmp;
+    THREAD_SAFE_STATIC multi_index mi;
     
     mi = t.sizes();
     size_type dimt = mi[ni], dim = m.ncols();
diff --git a/src/getfem/dal_basic.h b/src/getfem/dal_basic.h
index 878b996..932ca72 100644
--- a/src/getfem/dal_basic.h
+++ b/src/getfem/dal_basic.h
@@ -324,8 +324,8 @@ namespace dal
   template<class T, unsigned char pks> 
     typename dynamic_array<T,pks>::const_reference
       dynamic_array<T,pks>::operator [](size_type ii) const { 
-        DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(std::shared_ptr<T>,pf,NULL);
-        if (pf.get() == NULL) { pf = std::make_shared<T>(); }
+        THREAD_SAFE_STATIC std::shared_ptr<T> pf = nullptr;
+        if (pf.get() == nullptr) { pf = std::make_shared<T>(); }
         return (ii<last_ind) ? (array[ii>>pks])[ii&DNAMPKS__] : *pf;
   }
 
diff --git a/src/getfem/getfem_omp.h b/src/getfem/getfem_omp.h
index 17e8a7e..13d0691 100644
--- a/src/getfem/getfem_omp.h
+++ b/src/getfem/getfem_omp.h
@@ -47,8 +47,6 @@ This is the kernel of getfem.
 
 #ifdef GETFEM_HAS_OPENMP
   #include <mutex>
-  #include <boost/thread.hpp> /**TODO: get rid of this dependency as soon
-                                       as thread_local is widely supported*/
 #endif
 
 namespace getfem
@@ -341,33 +339,9 @@ namespace getfem
   and their initialization (it's more general and portable
   then using __declspec(thread))*/
   #ifdef GETFEM_HAS_OPENMP
-
-    #define DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(Type,Var,initial) \
-      static boost::thread_specific_ptr<Type> ptr_##Var; \
-      if(!ptr_##Var.get()) {ptr_##Var.reset(new Type(initial));} \
-      Type& Var=*ptr_##Var;
-
-    #define DEFINE_STATIC_THREAD_LOCAL(Type,Var) \
-      static boost::thread_specific_ptr<Type> ptr_##Var; \
-      if(!ptr_##Var.get()) {ptr_##Var.reset(new Type());} \
-      Type& Var=*ptr_##Var;
-
-    #define DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(Type, Var, ...) \
-      static boost::thread_specific_ptr<Type> ptr_##Var; \
-      if(!ptr_##Var.get()) {ptr_##Var.reset(new Type(__VA_ARGS__));} \
-      Type& Var=*ptr_##Var;
-
+    #define THREAD_SAFE_STATIC thread_local
   #else
-
-    #define DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(Type,Var,initial) \
-      static Type Var(initial);
-
-    #define DEFINE_STATIC_THREAD_LOCAL(Type,Var) \
-      static Type Var;
-
-    #define DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(Type, Var, ...) \
-      static Type Var(__VA_ARGS__);
-
+    #define THREAD_SAFE_STATIC static
   #endif
 
   class partition_master;
diff --git a/src/getfem_fem.cc b/src/getfem_fem.cc
index 263da5c..f0c0e9e 100644
--- a/src/getfem_fem.cc
+++ b/src/getfem_fem.cc
@@ -383,8 +383,8 @@ namespace getfem {
   typedef dal::dynamic_tree_sorted<dof_description, dof_description_comp__> 
dof_d_tab;
 
   pdof_description lagrange_dof(dim_type n) {
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(dim_type, n_old, dim_type(-2));
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pdof_description, p_old, 0);
+    THREAD_SAFE_STATIC dim_type n_old = dim_type(-2);
+    THREAD_SAFE_STATIC pdof_description p_old = nullptr;
     if (n != n_old) {
       dof_d_tab& tab = dal::singleton<dof_d_tab>::instance();
       dof_description l;
@@ -397,8 +397,8 @@ namespace getfem {
   }
 
   pdof_description lagrange_0_dof(dim_type n) {
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(dim_type, n_old, dim_type(-2));
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pdof_description, p_old, 0);
+    THREAD_SAFE_STATIC dim_type n_old = dim_type(-2);
+    THREAD_SAFE_STATIC pdof_description p_old = nullptr;
     if (n != n_old) {
       dof_d_tab& tab = dal::singleton<dof_d_tab>::instance();
       dof_description l;
@@ -3174,10 +3174,10 @@ namespace getfem {
   void hermite_segment__::mat_trans(base_matrix &M,
                                     const base_matrix &G,
                                     bgeot::pgeometric_trans pgt) const {
-    DEFINE_STATIC_THREAD_LOCAL(bgeot::pgeotrans_precomp, pgp);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(bgeot::pgeometric_trans, 
pgt_stored, 0);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_matrix, K, 1, 1);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_vector,r, 1);
+    THREAD_SAFE_STATIC bgeot::pgeotrans_precomp pgp;
+    THREAD_SAFE_STATIC bgeot::pgeometric_trans pgt_stored = nullptr;
+    THREAD_SAFE_STATIC base_matrix K(1, 1);
+    THREAD_SAFE_STATIC base_vector r(1);
     dim_type N = dim_type(G.nrows());
 
     if (pgt != pgt_stored) {
@@ -3240,9 +3240,9 @@ namespace getfem {
                                     const base_matrix &G,
                                     bgeot::pgeometric_trans pgt) const {
 
-    DEFINE_STATIC_THREAD_LOCAL(bgeot::pgeotrans_precomp, pgp);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(bgeot::pgeometric_trans, 
pgt_stored, 0);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_matrix, K, 2, 2);
+    THREAD_SAFE_STATIC bgeot::pgeotrans_precomp pgp;
+    THREAD_SAFE_STATIC bgeot::pgeometric_trans pgt_stored = nullptr;
+    THREAD_SAFE_STATIC base_matrix K(2, 2);
     dim_type N = dim_type(G.nrows());
 
     GMM_ASSERT1(N == 2, "Sorry, this version of hermite "
@@ -3312,9 +3312,9 @@ namespace getfem {
   void hermite_tetrahedron__::mat_trans(base_matrix &M,
                                     const base_matrix &G,
                                     bgeot::pgeometric_trans pgt) const {
-    DEFINE_STATIC_THREAD_LOCAL(bgeot::pgeotrans_precomp, pgp);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(bgeot::pgeometric_trans, 
pgt_stored, 0);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_matrix, K, 3, 3);
+    THREAD_SAFE_STATIC bgeot::pgeotrans_precomp pgp;
+    THREAD_SAFE_STATIC bgeot::pgeometric_trans pgt_stored = nullptr;
+    THREAD_SAFE_STATIC base_matrix K(3, 3);
     dim_type N = dim_type(G.nrows());
     GMM_ASSERT1(N == 3, "Sorry, this version of hermite "
                 "element works only on dimension three.")
@@ -3416,10 +3416,10 @@ namespace getfem {
                                     const base_matrix &G,
                                     bgeot::pgeometric_trans pgt) const {
 
-    DEFINE_STATIC_THREAD_LOCAL(bgeot::pgeotrans_precomp, pgp);
-    DEFINE_STATIC_THREAD_LOCAL(pfem_precomp, pfp);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(bgeot::pgeometric_trans, 
pgt_stored, 0);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_matrix, K, 2, 2);
+    THREAD_SAFE_STATIC bgeot::pgeotrans_precomp pgp;
+    THREAD_SAFE_STATIC pfem_precomp pfp;
+    THREAD_SAFE_STATIC bgeot::pgeometric_trans pgt_stored = nullptr;
+    THREAD_SAFE_STATIC base_matrix K(2, 2);
     dim_type N = dim_type(G.nrows());
     GMM_ASSERT1(N == 2, "Sorry, this version of argyris "
                 "element works only on dimension two.")
@@ -3458,7 +3458,7 @@ namespace getfem {
       M(5+6*k, 3+6*k) = c*c;     M(5+6*k, 4+6*k) = c*d;       M(5+6*k, 5+6*k) 
= d*d;
     }
 
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_matrix, W, 3, 21);
+    THREAD_SAFE_STATIC base_matrix W(3, 21);
     base_small_vector norient(M_PI, M_PI * M_PI);
     if (pgt->is_linear()) gmm::lu_inverse(K);
     for (unsigned i = 18; i < 21; ++i) {
@@ -3477,11 +3477,11 @@ namespace getfem {
       for (unsigned j = 0; j < 21; ++j)
         W(i-18, j) = t(j, 0, 0) * v[0] + t(j, 0, 1) * v[1];
     }
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_matrix,A,3,3);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(bgeot::base_vector, w, 3);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(bgeot::base_vector, coeff, 3);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(gmm::sub_interval, SUBI, 18,3);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(gmm::sub_interval, SUBJ, 0,3);
+    THREAD_SAFE_STATIC base_matrix A(3,3);
+    THREAD_SAFE_STATIC bgeot::base_vector w(3);
+    THREAD_SAFE_STATIC bgeot::base_vector coeff(3);
+    THREAD_SAFE_STATIC gmm::sub_interval SUBI(18, 3);
+    THREAD_SAFE_STATIC gmm::sub_interval SUBJ(0, 3);
     gmm::copy(gmm::sub_matrix(W, SUBJ, SUBI), A);
     gmm::lu_inverse(A);
     gmm::copy(gmm::transposed(A), gmm::sub_matrix(M, SUBI));
@@ -3587,10 +3587,10 @@ namespace getfem {
                                     const base_matrix &G,
                                     bgeot::pgeometric_trans pgt) const {
 
-    DEFINE_STATIC_THREAD_LOCAL(bgeot::pgeotrans_precomp, pgp);
-    DEFINE_STATIC_THREAD_LOCAL(pfem_precomp, pfp);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(bgeot::pgeometric_trans, 
pgt_stored, 0);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_matrix, K, 2, 2);
+    THREAD_SAFE_STATIC bgeot::pgeotrans_precomp pgp;
+    THREAD_SAFE_STATIC pfem_precomp pfp;
+    THREAD_SAFE_STATIC bgeot::pgeometric_trans pgt_stored = nullptr;
+    THREAD_SAFE_STATIC base_matrix K(2, 2);
     dim_type N = dim_type(G.nrows());
     GMM_ASSERT1(N == 2, "Sorry, this version of morley "
                 "element works only on dimension two.")
@@ -3601,7 +3601,7 @@ namespace getfem {
       pfp = fem_precomp(std::make_shared<morley_triangle__>(), node_tab(0), 0);
     }
     gmm::copy(gmm::identity_matrix(), M);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_matrix, W, 3, 6);
+    THREAD_SAFE_STATIC base_matrix W(3, 6);
     base_small_vector norient(M_PI, M_PI * M_PI);
     if (pgt->is_linear())
       { gmm::mult(G, pgp->grad(0), K); gmm::lu_inverse(K); }
@@ -3623,11 +3623,11 @@ namespace getfem {
     }
     //    cout << "W = " << W << endl; getchar();
 
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_matrix, A, 3, 3);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_vector, w, 3);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(base_vector, coeff, 3);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(gmm::sub_interval, SUBI, 3, 3);
-    DEFINE_STATIC_THREAD_LOCAL_CONSTRUCTED(gmm::sub_interval, SUBJ, 0, 3);
+    THREAD_SAFE_STATIC base_matrix A(3, 3);
+    THREAD_SAFE_STATIC base_vector w(3);
+    THREAD_SAFE_STATIC base_vector coeff(3);
+    THREAD_SAFE_STATIC gmm::sub_interval SUBI(3, 3);
+    THREAD_SAFE_STATIC gmm::sub_interval SUBJ(0, 3);
     gmm::copy(gmm::sub_matrix(W, SUBJ, SUBI), A);
     gmm::lu_inverse(A);
     gmm::copy(gmm::transposed(A), gmm::sub_matrix(M, SUBI));
@@ -3777,12 +3777,12 @@ namespace getfem {
   static pfem classical_fem_(bgeot::pgeometric_trans pgt,
                              short_type k, bool complete=false,
                              bool discont=false, scalar_type alpha=0) {
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(bgeot::pgeometric_trans, 
pgt_last,0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(short_type, k_last, short_type(-1));
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pfem, fem_last, 0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(char, complete_last, 0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(char, discont_last, 0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(scalar_type, alpha_last, 0);
+    THREAD_SAFE_STATIC bgeot::pgeometric_trans pgt_last = nullptr;
+    THREAD_SAFE_STATIC short_type k_last = short_type(-1);
+    THREAD_SAFE_STATIC pfem fem_last = nullptr;
+    THREAD_SAFE_STATIC char complete_last = 0;
+    THREAD_SAFE_STATIC char discont_last = 0;
+    THREAD_SAFE_STATIC scalar_type alpha_last = 0;
 
     bool found = false;
     if (pgt_last == pgt && k_last == k && complete_last == complete &&
@@ -3978,9 +3978,9 @@ namespace getfem {
   /* ******************************************************************** */
 
   pfem PK_fem(size_type n, short_type k) {
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pfem, pf, 0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(size_type, d, size_type(-2));
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(short_type, r, short_type(-2));
+    THREAD_SAFE_STATIC pfem pf = nullptr;
+    THREAD_SAFE_STATIC size_type d = size_type(-2);
+    THREAD_SAFE_STATIC short_type r = short_type(-2);
     if (d != n || r != k) {
       std::stringstream name;
       name << "FEM_PK(" << n << "," << k << ")";
@@ -3991,9 +3991,9 @@ namespace getfem {
   }
 
   pfem QK_fem(size_type n, short_type k) {
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pfem, pf, 0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(size_type, d, size_type(-2));
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(short_type, r, short_type(-2));
+    THREAD_SAFE_STATIC pfem pf = nullptr;
+    THREAD_SAFE_STATIC size_type d = size_type(-2);
+    THREAD_SAFE_STATIC short_type r = short_type(-2);
     if (d != n || r != k) {
       std::stringstream name;
       name << "FEM_QK(" << n << "," << k << ")";
@@ -4004,9 +4004,9 @@ namespace getfem {
   }
 
   pfem prism_PK_fem(size_type n, short_type k) {
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pfem, pf, 0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(size_type, d, size_type(-2));
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(short_type, r, short_type(-2));
+    THREAD_SAFE_STATIC pfem pf = nullptr;
+    THREAD_SAFE_STATIC size_type d = size_type(-2);
+    THREAD_SAFE_STATIC short_type r = short_type(-2);
     if (d != n || r != k) {
       std::stringstream name;
       name << "FEM_PRISM_PK(" << n << "," << k << ")";
diff --git a/src/getfem_integration.cc b/src/getfem_integration.cc
index 911a2b0..ce14aa5 100644
--- a/src/getfem_integration.cc
+++ b/src/getfem_integration.cc
@@ -1149,8 +1149,8 @@ namespace getfem {
   /* Fonctions pour la ref. directe.                                     */
 
   pintegration_method exact_simplex_im(size_type n) {
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pintegration_method, pim, 0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(size_type, d, -2);
+    THREAD_SAFE_STATIC pintegration_method pim = nullptr;
+    THREAD_SAFE_STATIC size_type d = -2;
     if (d != n) {
       std::stringstream name;
       name << "IM_EXACT_SIMPLEX(" << n << ")";
@@ -1161,8 +1161,8 @@ namespace getfem {
   }
 
   pintegration_method exact_parallelepiped_im(size_type n) {
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pintegration_method, pim, 0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(size_type, d, -2);
+    THREAD_SAFE_STATIC pintegration_method pim = nullptr;
+    THREAD_SAFE_STATIC size_type d = -2;
     if (d != n) {
       std::stringstream name;
       name << "IM_EXACT_PARALLELEPIPED(" << n << ")";
@@ -1173,8 +1173,8 @@ namespace getfem {
   }
 
   pintegration_method exact_prism_im(size_type n) {
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pintegration_method, pim, 0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(size_type, d, -2);
+    THREAD_SAFE_STATIC pintegration_method pim = nullptr;
+    THREAD_SAFE_STATIC size_type d = -2;
     if (d != n) {
       std::stringstream name;
       name << "IM_EXACT_PRISM(" << n << ")";
@@ -1190,8 +1190,8 @@ namespace getfem {
 
   static pintegration_method classical_exact_im(bgeot::pconvex_structure cvs) {
     cvs = bgeot::basic_structure(cvs);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(bgeot::pconvex_structure,cvs_last, 
0);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pintegration_method, im_last, 0);
+    THREAD_SAFE_STATIC bgeot::pconvex_structure cvs_last = nullptr;
+    THREAD_SAFE_STATIC pintegration_method im_last = nullptr;
     bool found = false;
 
     if (cvs_last == cvs)
@@ -1279,9 +1279,9 @@ namespace getfem {
 
   pintegration_method classical_approx_im(bgeot::pgeometric_trans pgt,
                                           dim_type degree) {
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(bgeot::pgeometric_trans, pgt_last, 
0);
-    DEFINE_STATIC_THREAD_LOCAL(dim_type, degree_last);
-    DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pintegration_method, im_last, 0);
+    THREAD_SAFE_STATIC bgeot::pgeometric_trans pgt_last = nullptr;
+    THREAD_SAFE_STATIC dim_type degree_last;
+    THREAD_SAFE_STATIC pintegration_method im_last = nullptr;
     if (pgt_last == pgt && degree == degree_last)
       return im_last;
     im_last = classical_approx_im_(pgt->structure(),degree);
@@ -1291,7 +1291,7 @@ namespace getfem {
   }
 
   pintegration_method im_none() {
-     DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pintegration_method,im_last,0);
+    THREAD_SAFE_STATIC pintegration_method im_last = nullptr;
      if (!im_last) im_last = int_method_descriptor("IM_NONE");
      return im_last;
   }
diff --git a/src/getfem_linearized_plates.cc b/src/getfem_linearized_plates.cc
index 55cba5e..d7fb4c3 100644
--- a/src/getfem_linearized_plates.cc
+++ b/src/getfem_linearized_plates.cc
@@ -105,8 +105,8 @@ namespace getfem {
     virtual void give_transformation(const mesh_fem &mf, size_type cv,
                                      base_matrix &M) const{
 
-      DEFINE_STATIC_THREAD_LOCAL(base_matrix, M_old);
-      DEFINE_STATIC_THREAD_LOCAL_INITIALIZED(pfem, pf_old, 0);
+      THREAD_SAFE_STATIC base_matrix M_old;
+      THREAD_SAFE_STATIC pfem pf_old = nullptr;
         
       // Obtaining the fem descriptors
       pfem pf1 = mf.fem_of_element(cv);
diff --git a/src/getfem_mat_elem.cc b/src/getfem_mat_elem.cc
index 9c5ea77..463552c 100644
--- a/src/getfem_mat_elem.cc
+++ b/src/getfem_mat_elem.cc
@@ -376,10 +376,10 @@ namespace getfem {
     void expand_product_daxpy(base_tensor &t, scalar_type J, bool first)const {
       size_type k;
       base_tensor::iterator pt = t.begin();
-      DEFINE_STATIC_THREAD_LOCAL(std::vector<base_tensor::const_iterator>, 
pts);
-      
DEFINE_STATIC_THREAD_LOCAL(std::vector<base_tensor::const_iterator>,es_beg);
-      
DEFINE_STATIC_THREAD_LOCAL(std::vector<base_tensor::const_iterator>,es_end);
-      DEFINE_STATIC_THREAD_LOCAL(std::vector<scalar_type>,Vtab);
+      THREAD_SAFE_STATIC std::vector<base_tensor::const_iterator> pts;
+      THREAD_SAFE_STATIC std::vector<base_tensor::const_iterator> es_beg;
+      THREAD_SAFE_STATIC std::vector<base_tensor::const_iterator> es_end;
+      THREAD_SAFE_STATIC std::vector<scalar_type> Vtab;
 
       pts.resize(0); pts.resize(pme->size()); // resize(0) necessary, do not 
remove
       es_beg.resize(0); es_beg.resize(pme->size());