diff --git a/CMake/FindFFTW.cmake b/CMake/FindFFTW.cmake
index 37f7c71527d8ee629e04b4141dda92f44007d84a..08ece4c55fcdbc1b5305b1d5770a97d333496f84 100644
--- a/CMake/FindFFTW.cmake
+++ b/CMake/FindFFTW.cmake
@@ -1,11 +1,63 @@
-# - Try to find FFTW
-# Once done, this will define
+# - Try to find FFTW3 libraries using components
 #
-#  FFTW_FOUND - system has FFTW
-#  FFTW_INCLUDE_DIRS - the FFTW include directories
-#  FFTW_LIBRARIES - link these to use FFTW
+# == Usage ==
+#   Call 
+#        find_package(FFTW 
+#             REQUIRED COMPONENTS <required components here>
+#            [OPTIONAL_COMPONENTS <optional components here>]
+#            [QUIET]
+#        )
+#   in your CmakeList.txt (note the missing underscore for required components)
+#   
+#   Possible component names: fftw3[fdlq](-(mpi|threads|omp))\?, where case doesn't matter
+#     and [fdlq] stands for f=float, d=double, l=long double, q=gcc quad long
+#   Note that the only exeption is fftw3q-mpi which doesn't exist at this day.
 #
-# Set fftw_DIR=where_fftw_is_installed if it's not in a "classic" place or if you want a specific version
+# == Full example == 
+#  set(FIND_FFTW_VERBOSE ON) (optional)
+#  set(FIND_FFTW_DEBUG OFF)  (optional)
+#  find_package(FFTW 
+#      REQUIRED COMPONENTS Fftw3f Fftw3d 
+#      OPTIONAL_COMPONENTS Fftw3l Fftw3q 
+#      OPTIONAL_COMPONENTS Fftw3f-mpi Fftw3d-mpi Fftw3l-mpi 
+#      OPTIONAL_COMPONENTS Fttw3f-omp Fftw3d-omp Fftw3l-omp Fftw3q-omp
+#      OPTIONAL_COMPONENTS Fttw3f-threads Fftw3d-threads Fftw3l-threads Fftw3q-threads)
+#
+# == Defined variables ==
+#
+#   FFTW_FOUND - will be set to true if and only if all required components and their dependencies were found
+#   FFTW_INCLUDE_DIRS - all FFTW include directories
+#   FFTW_LIBRARY_DIRS - all FFTW library directories
+#   FFTW_LIBRARIES - all FFTW libraries
+#   FFTW_DEFINES - all FFTW defines
+#   FTTW_COMPILE_FLAGS - all FFTW compile flags (currently only used for quadfloat fftw)
+#
+#   For each required or optional component, this package defines:
+#  		COMPONENT_FOUND   - will be set to true the component and its dependencies were found, else false
+#  		COMPONENT_DEFINES - all COMPONENT defines (-DFFTW_HAS_COMPONENT)
+#
+#  		COMPONENT_INCLUDE_DIRS - COMPONENT and its dependancies include directories
+#  		COMPONENT_LIBRARY_DIRS - COMPONENT and its dependancies library directories
+#  		COMPONENT_LIBRARIES    - COMPONENT and its dependancies libraries
+#
+#  		COMPONENT_INCLUDE_DIR - COMPONENT include directory
+#  		COMPONENT_LIBRARY_DIR - COMPONENT library directory
+#  		COMPONENT_LIBRARY     - COMPONENT library
+#
+#   where COMPONENT is equal to <input component name> transformed with the following rules:
+#       * uppercased
+#       * hyphens (-) replaced by underscores (_)
+#   Examples: fFtW3Q => FFTW3Q, fftw3f-mpi => FFTW3F_MPI
+# 
+# == Using a specific FFTW ==
+#   Set the variable ${fftw_DIR} to you desired search paths if it's not in a "classic" place or if you want a specific version. 
+#
+# == Checking against a specific version or the library ==
+#   Not supported yet.
+#
+# == Debug and verbose mode ==
+# Set ${FIND_FFTW_VERBOSE} to ON before find_package call to enable verbose mode
+# Set ${FIND_FFTW_DEBUG}   to ON before find_package call to enable debug mode
 #
 # Written by F. Pérignon, nov/2009
 # Updated by J-B. Keck, feb/2016
@@ -16,11 +68,6 @@ find_package(PkgConfig)
 # --version check not supported yet
 set(FFTW_VERSION_STRING ${FFTW_FIND_VERSION})
 
-# --main variables generated by this module
-list(APPEND FFTW_INCLUDE_DIRS "")
-list(APPEND FFTW_LIBRARY_DIRS "")
-list(APPEND FFTW_LIBRARIES "")
-        
 if(FIND_FFTW_VERBOSE)
     message(STATUS "Entering FindFFTW.cmake in verbose mode:")
 endif()
@@ -62,14 +109,16 @@ foreach(fftw_comp ${FFTW_FIND_COMPONENTS})
         NAMES ${header}
         PATHS ${fftw_DIR} 
         PATHS ${${COMPONENT}_PKGCONF_INCLUDE_DIRS}
-        PATHS ENV INCLUDE 
-              ENV PATH 
-              ENV C_INCLUDE_PATH 
         PATH_SUFFIXES include
         NO_DEFAULT_PATH
     )
     # -- search in default locations only if last search failed
-    find_path(FFTW_INCLUDE_DIR NAMES ${header})
+    find_path(${COMPONENT}_INCLUDE_DIR NAMES ${header}
+            PATHS ENV INCLUDE 
+                  ENV PATH 
+                  ENV C_INCLUDE_PATH 
+                  ENV CXX_INCLUDE_PATH 
+        )
 
     if(${${COMPONENT}_INCLUDE_DIR} STREQUAL "${COMPONENT}_INCLUDE_DIR-NOTFOUND")
         set(INCLUDE_DIR_FOUND FALSE)
@@ -90,14 +139,15 @@ foreach(fftw_comp ${FFTW_FIND_COMPONENTS})
       PATHS ${fftw_DIR} 
       PATHS ${${COMPONENT}_INCLUDE_DIR}/.. 
       PATHS ${${COMPONENT}_PKGCONF_LIBRARY_DIRS}}
-      PATHS ENV LIBRARY_PATH 
-            ENV LD_LIBRARY_PATH  
-            ENV DYLD_LIBRARY_PATH
       PATH_SUFFIXES lib
       NO_DEFAULT_PATH
     )
     # -- default locations
-    find_library(${COMPONENT}_LIBRARY NAMES ${library})
+    find_library(${COMPONENT}_LIBRARY 
+        NAMES ${library}
+        PATHS ENV LIBRARY_PATH 
+              ENV LD_LIBRARY_PATH  
+              ENV DYLD_LIBRARY_PATH)
     
     # -- if component is required append it to required vars
     if(FFTW_FIND_REQUIRED_${fftw_comp})
@@ -117,18 +167,29 @@ foreach(fftw_comp ${FFTW_FIND_COMPONENTS})
     # -- find quadmath library if required
     string(FIND ${component} "fftw3q" FFTWQ_POS)
     if(FFTWQ_POS EQUAL 0)
-        list(APPEND FFTW_LIBRARIES "quadmath")
+        if(NOT "${CMAKE_CXX_COMPILER_ID}" STREQUAL "GNU")
+            set(DEPENDENCIES_FOUND FALSE) # -- only gcc supports quadmaths
+        else()
+            list(APPEND ${COMPONENT}_LIBRARIES "quadmath")
+            list(APPEND ${COMPONENT}_DEFINES "-DHAS_QUADMATHS")  
+            list(APPEND FFTW_COMPILE_FLAGS "-fext-numeric-literals")
+            set(DEPENDENCIES_FOUND TRUE)
+        endif()
+    else()
+        set(DEPENDENCIES_FOUND TRUE)
     endif()
     
-    if(LIBRARY_DIR_FOUND AND INCLUDE_DIR_FOUND)
+    if(LIBRARY_DIR_FOUND AND INCLUDE_DIR_FOUND AND DEPENDENCIES_FOUND)
         set(${COMPONENT}_FOUND TRUE)
         list(APPEND ${COMPONENT}_INCLUDE_DIRS ${${COMPONENT}_INCLUDE_DIR})
         list(APPEND ${COMPONENT}_LIBRARY_DIRS ${${COMPONENT}_LIBRARY_DIR})
         list(APPEND ${COMPONENT}_LIBRARIES    ${${COMPONENT}_LIBRARY})
+        list(APPEND ${COMPONENT}_DEFINES "-DFFTW_HAS_${COMPONENT}")
 
         list(APPEND FFTW_INCLUDE_DIRS ${${COMPONENT}_INCLUDE_DIRS})
         list(APPEND FFTW_LIBRARY_DIRS ${${COMPONENT}_LIBRARY_DIRS})
         list(APPEND FFTW_LIBRARIES    ${${COMPONENT}_LIBRARIES})
+        list(APPEND FFTW_DEFINES      ${${COMPONENT}_DEFINES})
         
         if(FIND_FFTW_VERBOSE)
             message("\tFound FFTW::${fftw_comp} with parameters '-I${${COMPONENT}_INCLUDE_DIR} -L${${COMPONENT}_LIBRARY_DIR}  -l${${COMPONENT}_LIBRARY}'.")
@@ -149,9 +210,6 @@ foreach(fftw_comp ${FFTW_FIND_COMPONENTS})
         endif()
     endif()
         
-    unset(${COMPONENT}_INCLUDE_DIR)
-    unset(${COMPONENT}_LIBRARY_DIR)
-    unset(${COMPONENT}_LIBRARY)
     unset(FFTWQ_POS)
     unset(library)
     unset(header)
@@ -174,6 +232,7 @@ if(FIND_FFTW_DEBUG)
     message(STATUS "FFTW_INCLUDE_DIRS='${FFTW_INCLUDE_DIRS}'")
     message(STATUS "FFTW_LIBRARY_DIRS='${FFTW_LIBRARY_DIRS}'")
     message(STATUS "FFTW_LIBRARIES='${FFTW_LIBRARIES}'")
+    message(STATUS "FFTW_DEFINES='${FFTW_DEFINES}'")
 endif()
 
 unset(FFTW_REQUIRED_LIBRARIES)
diff --git a/CMake/HySoPInstallSetup.cmake b/CMake/HySoPInstallSetup.cmake
index 5a20120ca4084f04325b3029cc7e16754bb04126..ee14dcb52e25c91fcd9fbeb045cec30e954194d1 100644
--- a/CMake/HySoPInstallSetup.cmake
+++ b/CMake/HySoPInstallSetup.cmake
@@ -119,7 +119,3 @@ add_custom_target(python-install COMMAND ${PYTHON_EXECUTABLE} ${CMAKE_CURRENT_BI
 
 install(CODE "execute_process(COMMAND ${CMAKE_BUILD_TOOL} python-install WORKING_DIRECTORY \"${CMAKE_CURRENT_BINARY_DIR}\")")
 
-if(WITH_LIB_FORTRAN)
-  add_dependencies(python-install ${HYSOP_LIBRARY_NAME})
-
-endif()
diff --git a/HySoP/CMakeLists.txt b/HySoP/CMakeLists.txt
index 84253b2bdadb282a378405222495a8a7ab981693..fb0578d86cb96f9521d135cd42d970e964e2405b 100644
--- a/HySoP/CMakeLists.txt
+++ b/HySoP/CMakeLists.txt
@@ -55,7 +55,6 @@ set(hysop_python_install "user" CACHE STRING "Install mode for hysop python pack
 
 if(NOT WITH_LIB_FORTRAN)
   message(WARNING "You deactivate libhysop (fortran) generation. This will disable the fortran interface, including fftw and scales fonctionnalities.")
-  set(WITH_FFTW "OFF")
   set(WITH_SCALES "OFF")
   set(WITH_MAIN_FORTRAN "OFF")
 endif()
@@ -63,27 +62,31 @@ endif()
 if(NOT WITH_LIB_CXX)
     message(WARNING "You deactivate libhysop (cxx) generation. This will disable the Aitken-Schwarz Poisson solver.")
     set(WITH_MAIN_CXX "OFF")
+    set(WITH_GOOGLE_TESTS "OFF")
 endif()
 
-# Force a default build type if not provided by user
-# CMAKE_BUILD_TYPE = empty, Debug, Release, RelWithDebInfo or MinSizeRel.
-if(NOT CMAKE_BUILD_TYPE)
-  set(CMAKE_BUILD_TYPE RELEASE CACHE STRING
-    "Choose the type of build, options are: None Debug Release."
-    FORCE)
-endif()
-
-
 # true if hysop used Fortran and/or c++ sources
 if(${WITH_LIB_FORTRAN} OR ${WITH_LIB_CXX})
-  set(WITH_COMPILED_LIB)
+  set(WITH_COMPILED_LIB "ON")
+  set(WITH_FFTW "ON")
+else()
+  set(WITH_COMPILED_LIB "OFF")
+  set(WITH_FFTW "OFF")
 endif()
 
-# We can not run scales or fftw without mpi ...
+# We can not run scales without mpi ...
 if(WITH_FFTW OR WITH_SCALES)
   set(USE_MPI "ON")
 endif()
 
+# Force a default build type if not provided by user
+# CMAKE_BUILD_TYPE = empty, Debug, Release, RelWithDebInfo or MinSizeRel.
+if(NOT CMAKE_BUILD_TYPE)
+  set(CMAKE_BUILD_TYPE RELEASE CACHE STRING
+    "Choose the type of build, options are: None Debug Release."
+    FORCE)
+endif()
+
 # cmake project name
 set(PROJECT_NAME hysop)
 # --- Name for the package ---
@@ -119,14 +122,6 @@ set(PACKAGE_NAME HySoP)
 find_package(PythonFull REQUIRED)
 include(FindPythonModule)
 
-find_package(SWIG 3.0.8 REQUIRED)
-# WARNING FP : for cmake < 3.0 UseSWIG.cmake
-# does not work properly (bug for swig outdir)
-if(CMAKE_VERSION VERSION_LESS 3.0.0)
-    set(SWIG_USE_FILE ${CMAKE_SOURCE_DIR}/CMake/UseSWIG.cmake)
-endif()
-include(${SWIG_USE_FILE})
-
 # - python packages -
 find_python_module(numpy REQUIRED)
 find_python_module(scipy)
@@ -144,16 +139,6 @@ if(USE_MPI)
   find_python_module(mpi4py REQUIRED)
 endif()
 
-if(WITH_LIB_CXX)
-  find_package(SWIG 2.0.7 REQUIRED)
-  # WARNING FP : for cmake < 3.0 UseSWIG.cmake
-  # does not work properly (bug for swig outdir)
-  if(CMAKE_VERSION VERSION_LESS 3.0.0)
-    set(SWIG_USE_FILE ${CMAKE_SOURCE_DIR}/cmake/UseSWIG.cmake)
-  endif()
-  include(${SWIG_USE_FILE})
-endif()
-
 # Find python build dir name --> needed for tests and doc
 if(WITH_COMPILED_LIB)
   execute_process(
@@ -167,24 +152,24 @@ endif()
 
 # ============= Other dependencies =============
 
-# --- FFTW ---
-if(WITH_FFTW)
-  compile_with(FFTW REQUIRED)
-  set(dirlist)
-  foreach(_file ${FFTW_LIBRARIES})
-    if(CMAKE_PATCH_VERSION LESS 12)
-      get_filename_component(_name ${_file} PATH)
-    else()
-      get_filename_component(_name ${_file} DIRECTORY)
-    endif()
-    list(FIND dirlist ${_name} isin)
-    if(isin EQUAL -1)
-      list(APPEND dirlist ${_name})
+# -- Swig --
+if(WITH_LIB_CXX)
+    find_package(SWIG 3.0.8 REQUIRED)
+    if(CMAKE_VERSION VERSION_LESS 3.0.0)
+        set(SWIG_USE_FILE ${CMAKE_SOURCE_DIR}/CMake/UseSWIG.cmake)
     endif()
-  endforeach()
-  set(FFTWLIB ${dirlist} CACHE PATH "fftw libraries dir")
+    include(${SWIG_USE_FILE})
 endif()
 
+if(WITH_FFTW)
+    set(FIND_FFTW_VERBOSE OFF)
+    set(FIND_FFTW_DEBUG OFF)
+    find_package(FFTW 
+        REQUIRED COMPONENTS Fftw3d Fftw3f Fftw3d-mpi Fftw3f-mpi
+        OPTIONAL_COMPONENTS Fftw3l Fftw3q)
+endif()
+
+
 # ========= Check which opencl devices are available on the system =========
 if(WITH_GPU)
   execute_process(
@@ -248,27 +233,29 @@ if(WITH_LIB_FORTRAN)
 
   set(Fortran_FLAGS ${CMAKE_Fortran_FLAGS})
   append_flags(Fortran_FLAGS ${CMAKE_Fortran_FLAGS_${CMAKE_BUILD_TYPE}})
+endif()
 
+if(WITH_COMPILED_LIB)
   add_subdirectory(src)
+  #get_directory_property(FORTRAN_INCLUDE_DIRS
+    #DIRECTORY ${CMAKE_SOURCE_DIR}/src
+    #DEFINITION FORTRAN_INCLUDE_DIRS)
 endif()
 
 
 #C++ variables used by setup.py.in for swig
-set(CMAKE_CXX_FLAGS                "${CMAKE_CXX_FLAGS} -W -Wall -Wextra -Wno-unused-variable -Wno-unused-but-set-variable -Wno-unused-parameter -std=c++11")
+set(CMAKE_CXX_FLAGS                "${CMAKE_CXX_FLAGS} -W -Wall -Wextra -Wno-unused-variable -Wno-unused-but-set-variable -Wno-unused-parameter ${FFTW_COMPILE_FLAGS} -std=c++11")
 set(CMAKE_CXX_FLAGS_DEBUG          "${CMAKE_CXX_FLAGS_DEBUG}")
 set(CMAKE_CXX_FLAGS_RELWITHDEBINFO "${CMAKE_CXX_FLAGS_RELWITHDEBINFO}")
 set(CMAKE_CXX_FLAGS_RELEASE        "${CMAKE_CXX_FLAGS_RELEASE}")
 set(CMAKE_EXE_LINKER_FLAGS         "${CMAKE_EXE_LINKER_FLAGS}")
 
-if(NOT FFTW_FOUND)
-    find_package(FFTW REQUIRED)
-endif()
 set(CXX_FLAGS "${CMAKE_CXX_FLAGS}")
 set(CXX_LINKER_FLAGS "${CMAKE_EXE_LINKER_FLAGS}")
 set(CXX_EXT_INCLUDES ${PYTHON_INCLUDE_DIR} ${FFTW_INCLUDES})
 set(CXX_EXT_LIBS ${PYTHON_LIBRARIES} ${FFTW_LIBRARIES})
-set(CXX_EXT_LIB_DIRS "")
-set(CXX_EXTRA_DEFINES "")
+set(CXX_EXT_LIB_DIRS ${FFTW_LIBRARY_DIRS})
+set(CXX_EXTRA_DEFINES ${FFTW_DEFINES})
 
 #swig package name (lib name generated by swig)
 set(CPP_2_HYSOP "cpp2hysop")
@@ -294,8 +281,8 @@ set(PYTHON_SETUP "${CMAKE_CURRENT_BINARY_DIR}/setup.py")
 if(WITH_LIB_FORTRAN)
   add_custom_target(python-build ALL 
       COMMAND ${PYTHON_EXECUTABLE} ${PYTHON_SETUP} build config_fc --f90exec=${CMAKE_Fortran_COMPILER}
-      WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR} COMMENT "build hysop package")
-  add_dependencies(python-build ${HYSOP_LIBRARY_NAME})
+      WORKING_DIRECTORY ${CMAKE_CURRENT_BINARY_DIR} COMMENT "build hysop package"
+      DEPENDS ${HYSOP_LIBRARY_NAME}) #fortran module files have to be generated before
 else()
   add_custom_target(python-build ALL 
     COMMAND ${PYTHON_EXECUTABLE} ${PYTHON_SETUP} build
@@ -312,16 +299,16 @@ if(WITH_LIB_CXX)
         COMMAND cp `find ${CMAKE_CURRENT_BINARY_DIR}/build -name _${CPP_2_HYSOP}.so` ${HYSOP_CXX_LIBRARY}
         WORKING_DIRECTORY "${CMAKE_CURRENT_BINARY_DIR}/build"
         COMMENT "Copy swig c++ library to link")
-
+  
     get_filename_component(CXX_DIR  "${CMAKE_SOURCE_DIR}/src/hysop++" ABSOLUTE)
     get_filename_component(CXX_MAIN_DIR "${CXX_DIR}/main"             ABSOLUTE)
     get_filename_component(CXX_TEST_DIR "${CXX_DIR}/tests"            ABSOLUTE)
     get_filename_component(CXX_SOURCE_DIR "${CXX_DIR}/src"            ABSOLUTE)
     
-    list(APPEND cxx_local_include_dirs "${CXX_SOURCE_DIR}")
-    include_directories(${cxx_local_include_dirs})
+    include_directories(${CXX_SOURCE_DIR})
     include_directories(SYSTEM ${CXX_EXT_INCLUDES})
     link_directories(${CXX_EXT_LIB_DIRS})
+    add_definitions(${CXX_EXTRA_DEFINES})
 
     if(WITH_MAIN_CXX)
         list(APPEND cxx_executable_sources "${CXX_MAIN_DIR}/poissonSolver.cpp")
@@ -340,7 +327,6 @@ if(WITH_LIB_CXX)
 endif()
 
 
-
 # ====== Create a Target to clean sources (remove .pyc files) and build dir ======
 
 file(GLOB_RECURSE PYCFILES "${CMAKE_SOURCE_DIR}/*.pyc")
diff --git a/HySoP/setup.py.in b/HySoP/setup.py.in
index 64ddd44bf309adf9f0d1d53b265000c7ba254464..11d96582d35fc583c5d23e0073c6bffb655cc005 100644
--- a/HySoP/setup.py.in
+++ b/HySoP/setup.py.in
@@ -6,7 +6,9 @@ setup.py file for @PYPACKAGE_NAME@
 """
 from numpy.distutils.core import setup, Extension
 from numpy.distutils.misc_util import Configuration
+
 import os
+import re
 import glob
 import fnmatch
 
@@ -15,12 +17,30 @@ import sys
 sys.path.append('@CMAKE_SOURCE_DIR@/')
 import sort_f90
 
-def parseCmakeVar(var):
+def parseCMakeVar(var):
     if var != "":
         return var.split(';')
     else:
         return None
 
+def parseCMakeDefines(var):
+    defines = parseCMakeVar(var)
+    if defines == None:
+        return None
+   
+    # regex to match defines like -DMACRO_NAME = MACRO_VALUE
+    p = re.compile('\s*(?:-D)?\s*(\w+)(?:\s*=\s*(\w+))?\s*')
+
+    res = list()
+    for d in defines:
+        m = p.match(d)
+        if m:
+            res.append(m.group(1,2))
+        else:
+            print "\tWarning: Could extract cmake define from '",d,"'."
+    return res
+
+
 def create_fortran_extension(name, pyf_file=None, src_dirs=None, sources=None,
                              libdir=None, libs=None, debug_mode=0):
     """Create a new f2py module from fortran files
@@ -50,6 +70,12 @@ def create_fortran_extension(name, pyf_file=None, src_dirs=None, sources=None,
         inc_dir.remove('')
     inc_dir.append('@CMAKE_BINARY_DIR@/Modules')
     fortran_flags = ['@Fortran_FLAGS@']
+    
+    #includes = parseCMakeVar("@FORTRAN_INCLUDE_DIRS@")
+    #if(includes != None):
+        #for exti in includes:
+            #inc_dir.append(exti)
+
     ext_fort = Extension(name=name,
                          sources=sources,
                          f2py_options=f2py_options,
@@ -93,23 +119,23 @@ def create_swig_extension(name, inc_dirs, src_dirs=None, sources=None):
                  '-c++', '-modern']
                  #-outdir', '/Users/Franck/toto',
     
-    extern_includes = parseCmakeVar("@CXX_EXT_INCLUDES@")
+    extern_includes = parseCMakeVar("@CXX_EXT_INCLUDES@")
     if(extern_includes != None):
         for exti in extern_includes:
             include_dirs.append(exti)
 
-    libraries = parseCmakeVar("@CXX_EXT_LIBS@")
-    library_dirs = parseCmakeVar("@CXX_EXT_LIB_DIRS@")
-    extra_compile_args = parseCmakeVar("@CXX_FLAGS@")
-    extra_link_args = parseCmakeVar("@CXX_LINKER_FLAGS@")
-    define_macros = parseCmakeVar("@CXX_EXTRA_DEFINES@")
+    libraries = parseCMakeVar("@CXX_EXT_LIBS@")
+    library_dirs = parseCMakeVar("@CXX_EXT_LIB_DIRS@")
+    extra_compile_args = parseCMakeVar("@CXX_FLAGS@")
+    extra_link_args = parseCMakeVar("@CXX_LINKER_FLAGS@")
+    define_macros = parseCMakeDefines("@CXX_EXTRA_DEFINES@")
 
     #print "INCLUDE_DIRS=",include_dirs
     #print "LIBRARIES=",libraries
     #print "LIBRARY_DIRS=",library_dirs
     #print "EXTRA_COMPILE_ARGS=",extra_compile_args
     #print "EXTRA_LINK_ARGS=",extra_link_args
-    #print "DEFINE_MACROS=",define_macros
+    print "DEFINE_MACROS=",define_macros
     #  To avoid -I -I in compiler call, which results in a bug:
     #while inc_dir.count('') > 0:
     #    inc_dir.remove('')
diff --git a/HySoP/src/CMakeLists.txt b/HySoP/src/CMakeLists.txt
index 9364f9158499ed0a55a9c844c8de10fbc4a4c1ac..9d70b0d58211d8dd310d064c892632f8c7fa6c3f 100644
--- a/HySoP/src/CMakeLists.txt
+++ b/HySoP/src/CMakeLists.txt
@@ -8,23 +8,23 @@
 # The list of all dirs containing sources to be compiled for the HySoP lib
 # Any file in those dirs will be used to create libhysop
 set(${HYSOP_LIBRARY_NAME}_SRCDIRS
-  .
+    ${CMAKE_CURRENT_SOURCE_DIR}
   )
 
 # --- fftw interface ---
 if(WITH_FFTW)
   set(${HYSOP_LIBRARY_NAME}_SRCDIRS
-    ${${HYSOP_LIBRARY_NAME}_SRCDIRS} fftw
+    ${${HYSOP_LIBRARY_NAME}_SRCDIRS} ${CMAKE_CURRENT_SOURCE_DIR}/fftw
     )
 endif()
 if(WITH_LIB_CXX)
   set(${HYSOP_LIBRARY_NAME}_SRCDIRS
-    ${${HYSOP_LIBRARY_NAME}_SRCDIRS} hysop++/src
+    ${${HYSOP_LIBRARY_NAME}_SRCDIRS} ${CMAKE_CURRENT_SOURCE_DIR}/hysop++/src
     )
 endif()
 
 #set(SCALES_DIR scalesReduced)
-set(SCALES_DIR scalesInterface)
+set(SCALES_DIR ${CMAKE_CURRENT_SOURCE_DIR}/scalesInterface)
 
 # --- scales ---
 if(WITH_SCALES)
@@ -36,6 +36,7 @@ if(WITH_SCALES)
     ${SCALES_DIR}/particles/advec_line
     )
 endif()
+list(APPEND FORTRAN_INCLUDE_DIRS ${${HYSOP_LIBRARY_NAME}_SRCDIRS})
 
 # Matching expr for files to be compiled.
 set(EXTS *.f90 *.f95 *F90)
@@ -69,6 +70,7 @@ set(CMAKE_Fortran_MODULE_DIRECTORY ${CMAKE_BINARY_DIR}/Modules)
 #append_Fortran_FLAGS("-Wall -fPIC -ffree-line-length-none -DBLOCKING_SEND_PLUS -DBLOCKING_SEND")
 append_Fortran_FLAGS("-Wall -fPIC -ffree-line-length-none")
 
+# --- MPI ---
 if(USE_MPI)
   # Find MPI for fortran.
   find_package(MPI REQUIRED)
@@ -83,26 +85,13 @@ endif(USE_MPI)
 
 # --- FFTW ---
 if(WITH_FFTW)
-  find_package(FFTW REQUIRED)
-  include_directories(${FFTW_INCLUDE_DIRS})
+  include_directories(SYSTEM ${FFTW_INCLUDE_DIRS})
+  link_directories(${FFTW_LIBRARY_DIRS})
+  add_definitions(${FFTW_DEFINES})
   set(LIBS ${LIBS} ${FFTW_LIBRARIES})
-  # for python setup.py
-  #set(FFTWLIB ${FFTW_LIBRARIES} CACHE PATH "fftw libraries dir")
-  set(dirlist)
-  foreach(_file ${FFTW_LIBRARIES})
-    if(CMAKE_PATCH_VERSION LESS 12)
-      get_filename_component(_name ${_file} PATH)
-    else()
-      get_filename_component(_name ${_file} DIRECTORY)
-    endif()
-    list(FIND dirlist ${_name} isin)
-    if(isin EQUAL -1)
-      list(APPEND dirlist ${_name})
-    endif()
-  endforeach()
-  set(FFTWLIB ${dirlist} CACHE PATH "fftw libraries dir")
 endif()
 
+
 # ============= Generates HySoPConfig.hpp =============
 # The file HYSOP_LIBRARY_NAME_defines.hpp will be generated from config.hpp.cmake;
 if(EXISTS ${CMAKE_SOURCE_DIR}/config.hpp.cmake)
@@ -132,8 +121,7 @@ list(APPEND ${HYSOP_LIBRARY_NAME}_HDRS "hysop_defines.hpp")
 
 # We add headers to source files list
 list(APPEND ${HYSOP_LIBRARY_NAME}_SRC ${${HYSOP_LIBRARY_NAME}_HDRS})
-# Add directories to those searched by compiler ...
-# -I
+
 include_directories(${${HYSOP_LIBRARY_NAME}_SRCDIRS})
 include_directories(${CMAKE_Fortran_MODULE_DIRECTORY})
 
diff --git a/HySoP/src/hysop++/main/diffSolver.cpp b/HySoP/src/hysop++/main/diffSolver.cpp
index 0acad64562fd65552a5335d25c44682177121a98..5cf98f3c9c4df4200ced9e09888325ac54f53b81 100644
--- a/HySoP/src/hysop++/main/diffSolver.cpp
+++ b/HySoP/src/hysop++/main/diffSolver.cpp
@@ -6,6 +6,7 @@
 #include "data/multi_array/multi_array.h"
 #include "utils/constants.h"
 #include "fft/extension.h"
+#include "maths/quad_maths.h"
 
 using namespace hysop;
 using namespace hysop::domain;
@@ -188,6 +189,7 @@ void test(std::size_t p_maxOrder, bool includePeriodicBds=false) {
     
 int main(int argc, const char *argv[]) {
 
+#ifdef FFTW_HAS_FFTW3F
     std::cout << "== TEST 1D - float       ==" << std::endl;
     test<float,1,false>(5);
     std::cout << "== TEST 2D - float       ==" << std::endl;
@@ -195,7 +197,9 @@ int main(int argc, const char *argv[]) {
     std::cout << "== TEST 3D - float       ==" << std::endl;
     test<float,3,false>(1);
     std::cout << std::endl;
+#endif
     
+#ifdef FFTW_HAS_FFTW3D
     std::cout << "== TEST 1D - double       ==" << std::endl;
     test<double,1,false>(5);
     std::cout << "== TEST 2D - double       ==" << std::endl;
@@ -203,7 +207,9 @@ int main(int argc, const char *argv[]) {
     std::cout << "== TEST 3D - double       ==" << std::endl;
     test<double,3,false>(1);
     std::cout << std::endl;
+#endif
     
+#ifdef FFTW_HAS_FFTW3L
     std::cout << "== TEST 1D - long double       ==" << std::endl;
     test<long double,1,false>(5);
     std::cout << "== TEST 2D - long double       ==" << std::endl;
@@ -211,6 +217,17 @@ int main(int argc, const char *argv[]) {
     std::cout << "== TEST 3D - long double       ==" << std::endl;
     test<long double,3,false>(1);
     std::cout << std::endl;
+#endif
+
+//#ifdef FFTW_HAS_FFTW3Q
+    //std::cout << "== TEST 1D - __float128       ==" << std::endl;
+    //test<__float128,1,false>(5);
+    //std::cout << "== TEST 2D - __float128       ==" << std::endl;
+    //test<__float128,2,false>(3);
+    //std::cout << "== TEST 3D - __float128       ==" << std::endl;
+    //test<__float128,3,false>(1);
+    //std::cout << std::endl;
+//#endif
 
     return EXIT_SUCCESS;
 }
diff --git a/HySoP/src/hysop++/src/data/multi_array/multi_array.h b/HySoP/src/hysop++/src/data/multi_array/multi_array.h
index 6224ac4ca54709aac895138e80c4356292a197ec..41746a7eb7a451708e6aaaae88b297b1a4d6a099 100644
--- a/HySoP/src/hysop++/src/data/multi_array/multi_array.h
+++ b/HySoP/src/hysop++/src/data/multi_array/multi_array.h
@@ -2,8 +2,9 @@
 #ifndef HYSOP_MULTI_ARRAY_H
 #define HYSOP_MULTI_ARRAY_H
 
-#include <boost/multi_array.hpp>
+#include "utils/utils.h"
 #include "utils/default.h"
+#include <boost/multi_array.hpp>
 
 /****************************************/
 /*** Hysop boost multi array wrapper ****/
diff --git a/HySoP/src/hysop++/src/data/multi_array/multi_array_defines.h b/HySoP/src/hysop++/src/data/multi_array/multi_array_defines.h
index c28bd8ec013d7423a212884cc90ea7d224a965a3..75e51a227a59164f71db11addd7e36e238b39118 100644
--- a/HySoP/src/hysop++/src/data/multi_array/multi_array_defines.h
+++ b/HySoP/src/hysop++/src/data/multi_array/multi_array_defines.h
@@ -45,23 +45,15 @@
             typename std::enable_if<                                     \
                        std::is_same<TT,bool>::value                      \
                 , TYPE>::type
-#define ENABLE_IF_REAL(TYPE,DEFAULT)                                     \
+#define ENABLE_IF_FFTW_REAL(TYPE,DEFAULT)                                \
             template <typename TT DEFAULT>                               \
             typename std::enable_if<                                     \
-                       std::is_same<TT,float>      ::value               \
-                    || std::is_same<TT,double>     ::value               \
-                    || std::is_same<TT,long double>::value               \
-                , TYPE>::type
-#define ENABLE_IF_COMPLEX(TYPE,DEFAULT)                                  \
-            template <typename TT DEFAULT>                               \
-            typename std::enable_if<                                     \
-                       std::is_same<TT,std::complex<float>>      ::value \
-                    || std::is_same<TT,std::complex<double>>     ::value \
-                    || std::is_same<TT,std::complex<long double>>::value \
-                    || std::is_same<TT,fftwf_complex>            ::value \
-                    || std::is_same<TT, fftw_complex>            ::value \
-                    || std::is_same<TT,fftwl_complex>            ::value \
-                    || std::is_same<TT,fftwq_complex>            ::value \
+                    hysop::fft::is_fftw_supported_type<TT>::value        \
+                , TYPE>::type 
+#define ENABLE_IF_FFTW_COMPLEX(TYPE,DEFAULT)                               \
+            template <typename TT DEFAULT>                                 \
+            typename std::enable_if<                                       \
+                    hysop::fft::is_fftw_supported_complex_type<TT>::value  \
                 , TYPE>::type
 
 /* CLASS INTERFACES */
@@ -85,13 +77,13 @@
         /* common const interface */                                                                            \
         PUBLIC_COMMON_CONST_INTERFACE(SINGLE_ARG(TYPENAME))                                                      \
         /* real and complex data accessors, usefull for FFT like transforms */                                  \
-        ENABLE_IF_REAL   (const T*,=T) rdata() const;                                                           \
-        ENABLE_IF_REAL   (const typename fft::fftw_complex_type<T>::std_type*,=T)  asStdComplexData() const;    \
-        ENABLE_IF_REAL   (const typename fft::fftw_complex_type<T>::fftw_type*,=T) asFftwComplexData() const;   \
-        ENABLE_IF_COMPLEX(const T*,=T) cdata() const;                                                           \
-        ENABLE_IF_COMPLEX(const typename fft::fftw_complex_type<T>::std_type*,=T)    std_cdata() const;         \
-        ENABLE_IF_COMPLEX(const typename fft::fftw_complex_type<T>::fftw_type*,=T)  fftw_cdata() const;         \
-        ENABLE_IF_COMPLEX(const typename fft::fftw_complex_type<T>::value_type*,=T) asRealData() const;
+        ENABLE_IF_FFTW_REAL   (const T*,=T) rdata() const;                                                           \
+        ENABLE_IF_FFTW_REAL   (const typename fft::fftw_complex_type<T>::std_type*,=T)  asStdComplexData() const;    \
+        ENABLE_IF_FFTW_REAL   (const typename fft::fftw_complex_type<T>::fftw_type*,=T) asFftwComplexData() const;   \
+        ENABLE_IF_FFTW_COMPLEX(const T*,=T) cdata() const;                                                           \
+        ENABLE_IF_FFTW_COMPLEX(const typename fft::fftw_complex_type<T>::std_type*,=T)    std_cdata() const;         \
+        ENABLE_IF_FFTW_COMPLEX(const typename fft::fftw_complex_type<T>::fftw_type*,=T)  fftw_cdata() const;         \
+        ENABLE_IF_FFTW_COMPLEX(const typename fft::fftw_complex_type<T>::value_type*,=T) asRealData() const;
 
 
 /* public const view interface */
@@ -113,13 +105,13 @@
         /* common non const interface */                                                          \
         PUBLIC_COMMON_NON_CONST_INTERFACE(SINGLE_ARG(TYPENAME))                                   \
         /* real and complex data accessors, usefull for FFT like transforms */                    \
-        ENABLE_IF_REAL(T*,=T) rdata();                                                            \
-        ENABLE_IF_REAL(typename fft::fftw_complex_type<T>::std_type*,=T) asStdComplexData();      \
-        ENABLE_IF_REAL(typename fft::fftw_complex_type<T>::fftw_type*,=T) asFftwComplexData();    \
-        ENABLE_IF_COMPLEX(T*,=T) cdata();                                                         \
-        ENABLE_IF_COMPLEX(typename fft::fftw_complex_type<T>::std_type*,=T) std_cdata();          \
-        ENABLE_IF_COMPLEX(typename fft::fftw_complex_type<T>::fftw_type*,=T) fftw_cdata();        \
-        ENABLE_IF_COMPLEX(typename fft::fftw_complex_type<T>::value_type*,=T) asRealData();
+        ENABLE_IF_FFTW_REAL(T*,=T) rdata();                                                            \
+        ENABLE_IF_FFTW_REAL(typename fft::fftw_complex_type<T>::std_type*,=T) asStdComplexData();      \
+        ENABLE_IF_FFTW_REAL(typename fft::fftw_complex_type<T>::fftw_type*,=T) asFftwComplexData();    \
+        ENABLE_IF_FFTW_COMPLEX(T*,=T) cdata();                                                         \
+        ENABLE_IF_FFTW_COMPLEX(typename fft::fftw_complex_type<T>::std_type*,=T) std_cdata();          \
+        ENABLE_IF_FFTW_COMPLEX(typename fft::fftw_complex_type<T>::fftw_type*,=T) fftw_cdata();        \
+        ENABLE_IF_FFTW_COMPLEX(typename fft::fftw_complex_type<T>::value_type*,=T) asRealData();
 
 
 #define PUBLIC_NON_CONST_VIEW_INTERFACE(TYPENAME)                                                 \
@@ -246,31 +238,31 @@
 /* Reference specific const implementation */
 #define CONST_REF_IMPL(TYPENAME,TEMPLATES)                                                                                          \
     TEMPLATES                                                                                                                       \
-    ENABLE_IF_REAL(const T*,NO_DEFAULT_TEMPLATES) TYPENAME::rdata() const {                                                         \
+    ENABLE_IF_FFTW_REAL(const T*,NO_DEFAULT_TEMPLATES) TYPENAME::rdata() const {                                                         \
         return this->data();                                                                                                        \
     }                                                                                                                               \
     TEMPLATES                                                                                                                       \
-    ENABLE_IF_REAL(const typename fft::fftw_complex_type<T>::std_type*,NO_DEFAULT_TEMPLATES)  TYPENAME::asStdComplexData() const {  \
+    ENABLE_IF_FFTW_REAL(const typename fft::fftw_complex_type<T>::std_type*,NO_DEFAULT_TEMPLATES)  TYPENAME::asStdComplexData() const {  \
         return reinterpret_cast<const typename fft::fftw_complex_type<T>::std_type*>(this->data());                                 \
     }                                                                                                                               \
     TEMPLATES                                                                                                                       \
-    ENABLE_IF_REAL(const typename fft::fftw_complex_type<T>::fftw_type*,NO_DEFAULT_TEMPLATES) TYPENAME::asFftwComplexData() const { \
+    ENABLE_IF_FFTW_REAL(const typename fft::fftw_complex_type<T>::fftw_type*,NO_DEFAULT_TEMPLATES) TYPENAME::asFftwComplexData() const { \
         return reinterpret_cast<const typename fft::fftw_complex_type<T>::fftw_type*>(this->data());                                \
     }                                                                                                                               \
     TEMPLATES                                                                                                                       \
-    ENABLE_IF_COMPLEX(const T*,NO_DEFAULT_TEMPLATES) TYPENAME::cdata() const {                                                      \
+    ENABLE_IF_FFTW_COMPLEX(const T*,NO_DEFAULT_TEMPLATES) TYPENAME::cdata() const {                                                      \
         return this->data();                                                                                                        \
     }                                                                                                                               \
     TEMPLATES                                                                                                                       \
-    ENABLE_IF_COMPLEX(const typename fft::fftw_complex_type<T>::std_type*,NO_DEFAULT_TEMPLATES)    TYPENAME::std_cdata() const {    \
+    ENABLE_IF_FFTW_COMPLEX(const typename fft::fftw_complex_type<T>::std_type*,NO_DEFAULT_TEMPLATES)    TYPENAME::std_cdata() const {    \
         return reinterpret_cast<const typename fft::fftw_complex_type<T>::std_type*>(this->data());                                 \
     }                                                                                                                               \
     TEMPLATES                                                                                                                       \
-    ENABLE_IF_COMPLEX(const typename fft::fftw_complex_type<T>::fftw_type*,NO_DEFAULT_TEMPLATES)  TYPENAME::fftw_cdata() const {    \
+    ENABLE_IF_FFTW_COMPLEX(const typename fft::fftw_complex_type<T>::fftw_type*,NO_DEFAULT_TEMPLATES)  TYPENAME::fftw_cdata() const {    \
         return reinterpret_cast<const typename fft::fftw_complex_type<T>::fftw_type*>(this->data());                                \
     }                                                                                                                               \
     TEMPLATES                                                                                                                       \
-    ENABLE_IF_COMPLEX(const typename fft::fftw_complex_type<T>::value_type*,NO_DEFAULT_TEMPLATES) TYPENAME::asRealData() const {    \
+    ENABLE_IF_FFTW_COMPLEX(const typename fft::fftw_complex_type<T>::value_type*,NO_DEFAULT_TEMPLATES) TYPENAME::asRealData() const {    \
         return reinterpret_cast<const typename fft::fftw_complex_type<T>::value_type*>(this->data());                               \
     }                                                                                                                               
 
@@ -325,31 +317,31 @@
 /* Reference specific non const implementation */
 #define NON_CONST_REF_IMPL(TYPENAME,TEMPLATES)                                                                          \
     TEMPLATES                                                                                                           \
-    ENABLE_IF_REAL(T*,NO_DEFAULT_TEMPLATES) TYPENAME::rdata() {                                                         \
+    ENABLE_IF_FFTW_REAL(T*,NO_DEFAULT_TEMPLATES) TYPENAME::rdata() {                                                         \
         return this->data();                                                                                            \
     }                                                                                                                   \
     TEMPLATES                                                                                                           \
-    ENABLE_IF_REAL(typename fft::fftw_complex_type<T>::std_type*,NO_DEFAULT_TEMPLATES)  TYPENAME::asStdComplexData() {  \
+    ENABLE_IF_FFTW_REAL(typename fft::fftw_complex_type<T>::std_type*,NO_DEFAULT_TEMPLATES)  TYPENAME::asStdComplexData() {  \
         return reinterpret_cast<typename fft::fftw_complex_type<T>::std_type*>(this->data());                           \
     }                                                                                                                   \
     TEMPLATES                                                                                                           \
-    ENABLE_IF_REAL(typename fft::fftw_complex_type<T>::fftw_type*,NO_DEFAULT_TEMPLATES) TYPENAME::asFftwComplexData() { \
+    ENABLE_IF_FFTW_REAL(typename fft::fftw_complex_type<T>::fftw_type*,NO_DEFAULT_TEMPLATES) TYPENAME::asFftwComplexData() { \
         return reinterpret_cast<typename fft::fftw_complex_type<T>::fftw_type*>(this->data());                          \
     }                                                                                                                   \
     TEMPLATES                                                                                                           \
-    ENABLE_IF_COMPLEX(T*,NO_DEFAULT_TEMPLATES) TYPENAME::cdata() {                                                      \
+    ENABLE_IF_FFTW_COMPLEX(T*,NO_DEFAULT_TEMPLATES) TYPENAME::cdata() {                                                      \
         return this->data();                                                                                            \
     }                                                                                                                   \
     TEMPLATES                                                                                                           \
-    ENABLE_IF_COMPLEX(typename fft::fftw_complex_type<T>::std_type*,NO_DEFAULT_TEMPLATES)    TYPENAME::std_cdata() {    \
+    ENABLE_IF_FFTW_COMPLEX(typename fft::fftw_complex_type<T>::std_type*,NO_DEFAULT_TEMPLATES)    TYPENAME::std_cdata() {    \
         return reinterpret_cast<typename fft::fftw_complex_type<T>::std_type*>(this->data());                           \
     }                                                                                                                   \
     TEMPLATES                                                                                                           \
-    ENABLE_IF_COMPLEX(typename fft::fftw_complex_type<T>::fftw_type*,NO_DEFAULT_TEMPLATES)  TYPENAME::fftw_cdata() {    \
+    ENABLE_IF_FFTW_COMPLEX(typename fft::fftw_complex_type<T>::fftw_type*,NO_DEFAULT_TEMPLATES)  TYPENAME::fftw_cdata() {    \
         return reinterpret_cast<typename fft::fftw_complex_type<T>::fftw_type*>(this->data());                          \
     }                                                                                                                   \
     TEMPLATES                                                                                                           \
-    ENABLE_IF_COMPLEX(typename fft::fftw_complex_type<T>::value_type*,NO_DEFAULT_TEMPLATES) TYPENAME::asRealData() {    \
+    ENABLE_IF_FFTW_COMPLEX(typename fft::fftw_complex_type<T>::value_type*,NO_DEFAULT_TEMPLATES) TYPENAME::asRealData() {    \
         return reinterpret_cast<typename fft::fftw_complex_type<T>::value_type*>(this->data());                         \
     }                                                                                                                               
 
diff --git a/HySoP/src/hysop++/src/fft/fftw3.h b/HySoP/src/hysop++/src/fft/fftw3.h
index 6d214c5a4f9603b98af301f054dbc02c49b5e0fd..f08bca38b278a6e1a7417dee1c27e4082b867493 100644
--- a/HySoP/src/hysop++/src/fft/fftw3.h
+++ b/HySoP/src/hysop++/src/fft/fftw3.h
@@ -5,8 +5,8 @@
 #include <complex>
 #include <fftw3.h>
 
-#ifdef HYSOP_ENABLE_QUAD_MATH
-#include <quadmath.h>
+#ifdef HAS_QUADMATHS
+#include "maths/quad_maths.h"
 #endif
 
 /*                                                         */
@@ -403,30 +403,34 @@ namespace hysop {
         template <> struct is_fftw_supported_type<float>       { static constexpr bool value = true;  };
         template <> struct is_fftw_supported_type<double>      { static constexpr bool value = true;  };
         template <> struct is_fftw_supported_type<long double> { static constexpr bool value = true;  };
-#ifdef HYSOP_ENABLE_QUAD_MATH
-        template <> struct is_fftw_supported_type<__float128>  { static constexpr bool value = true;  };
-#endif
+        
+        template <typename T> struct is_fftw_supported_complex_type    { static constexpr bool value = false; };
+        template <> struct is_fftw_supported_complex_type<std::complex<float>>       { static constexpr bool value = true;  };
+        template <> struct is_fftw_supported_complex_type<std::complex<double>>      { static constexpr bool value = true;  };
+        template <> struct is_fftw_supported_complex_type<std::complex<long double>> { static constexpr bool value = true;  };
+        template <> struct is_fftw_supported_complex_type<fftwf_complex>             { static constexpr bool value = true;  };
+        template <> struct is_fftw_supported_complex_type<fftw_complex>              { static constexpr bool value = true;  };
+        template <> struct is_fftw_supported_complex_type<fftwl_complex>             { static constexpr bool value = true;  };
 
         template <typename T>
             struct Fftw3 {
                 Fftw3() { 
                     throw std::runtime_error(
                             "Can only use Fftw3 wrapper with types {float, double, long double, __float128} ! "
-                            "Special __float128 type is enabled via the macro HYSOP_ENABLE_QUAD_MATH."
+                            "Note: __float128 type is enabled only if HAS_QUADMATHS is defined."
                             );
                 }
             };
        
         /* Generate Ffftw<> template specialisations */
-        /* (1) link with lib fftw3f (-lfftw3f) */
-        /* (2) link with lib fftw3  (-lfftw3)  */
-        /* (3) link with lib fftw3l (-lfftw3l) */
-        /* (4) link with lib fftw3q and quadmath (-lfftw3q -lquadmath) */
-
         FFTW_DEFINE_CXX_API(FFTW_MANGLE_CLASS, FFTW_MANGLE_FLOAT, float)
         FFTW_DEFINE_CXX_API(FFTW_MANGLE_CLASS, FFTW_MANGLE_DOUBLE, double)
         FFTW_DEFINE_CXX_API(FFTW_MANGLE_CLASS, FFTW_MANGLE_LONG_DOUBLE, long double)
-#ifdef HYSOP_ENABLE_QUAD_MATH
+
+#ifdef HAS_QUADMATHS
+        template <> struct is_fftw_supported_type<__float128>                        { static constexpr bool value = true;  };
+        template <> struct is_fftw_supported_complex_type<std::complex<__float128>>  { static constexpr bool value = true;  };
+        template <> struct is_fftw_supported_complex_type<fftwq_complex>             { static constexpr bool value = true;  };
         FFTW_DEFINE_CXX_API(FFTW_MANGLE_CLASS, FFTW_MANGLE_QUAD, __float128)
 #endif
         
diff --git a/HySoP/src/hysop++/src/fft/fftwComplex.h b/HySoP/src/hysop++/src/fft/fftwComplex.h
index 711be112d779d5a0da9124b6cf348d5ccd116016..e2fd44593400408297fc02f397575e16b753eab7 100644
--- a/HySoP/src/hysop++/src/fft/fftwComplex.h
+++ b/HySoP/src/hysop++/src/fft/fftwComplex.h
@@ -4,16 +4,12 @@
 
 #include <type_traits>
 #include <complex>
-#include <fftw3.h>
+#include <fft/fftw3.h>
 
 namespace hysop {
     namespace fft {
 
-        template <typename T> struct fftw_complex_type { 
-            typedef void value_type;
-            typedef void fftw_type;
-            typedef void std_type;
-        };
+        template <typename T> struct fftw_complex_type {};
 
         template <> struct fftw_complex_type<long double> { 
             typedef long double               value_type;
@@ -63,6 +59,24 @@ namespace hysop {
             typedef std::complex<float> std_type; 
         };
 
+#ifdef HAS_QUADMATHS
+        template <> struct fftw_complex_type<__float128> { 
+            typedef __float128               value_type;
+            typedef fftwq_complex            fftw_type; 
+            typedef std::complex<__float128> std_type; 
+        };
+        template <> struct fftw_complex_type<std::complex<__float128>> { 
+            typedef __float128               value_type;
+            typedef fftwq_complex            fftw_type; 
+            typedef std::complex<__float128> std_type; 
+        };
+        template <> struct fftw_complex_type<fftwq_complex> { 
+            typedef __float128               value_type;
+            typedef fftwq_complex            fftw_type; 
+            typedef std::complex<__float128> std_type; 
+        };
+#endif
+
     }
 }
 
diff --git a/HySoP/src/hysop++/src/maths/quad_maths.h b/HySoP/src/hysop++/src/maths/quad_maths.h
new file mode 100644
index 0000000000000000000000000000000000000000..d02b3da6071af710ea9444c1fa9d892c28ce6712
--- /dev/null
+++ b/HySoP/src/hysop++/src/maths/quad_maths.h
@@ -0,0 +1,291 @@
+
+#ifdef HAS_QUADMATHS
+
+#ifndef HYSOP_QUAD_MATHS_H
+#define HYSOP_QUAD_MATHS_H
+
+#include <cfloat>
+#include <limits>
+#include <iostream>
+#include <iomanip>
+
+#include <quadmath.h>
+
+/* missing gcc defines */
+#define FLT128_RADIX FLT_RADIX
+#define FLT128_HAS_DENORM      true
+#define FLT128_HAS_INFINITY    true
+#define FLT128_HAS_QUIET_NAN   true
+
+namespace std {
+
+    template<> struct numeric_limits<__float128> {
+        static constexpr bool is_specialized = true;
+
+        static constexpr __float128 min() { return FLT128_MIN; }
+        static constexpr __float128 max() { return FLT128_MAX; }
+        static constexpr __float128 lowest() noexcept { return -FLT128_MAX; }
+
+        static constexpr int digits = FLT128_MANT_DIG;
+        static constexpr int digits10 = FLT128_DIG;
+        static constexpr int max_digits10 = (2 + (FLT128_MANT_DIG) * 643L / 2136);
+
+        static constexpr bool is_signed = true;
+        static constexpr bool is_integer = false;
+        static constexpr bool is_exact = false;
+        static constexpr int radix = FLT128_RADIX;
+
+        static constexpr __float128 epsilon()      { return FLT128_EPSILON; }
+        static constexpr __float128 round_error()  { return 0.5; }
+
+        static constexpr int min_exponent = FLT128_MIN_EXP;
+        static constexpr int min_exponent10 = FLT128_MIN_10_EXP;
+        static constexpr int max_exponent = FLT128_MAX_EXP;
+        static constexpr int max_exponent10 = FLT128_MAX_10_EXP;
+
+        static constexpr bool has_infinity = FLT128_HAS_INFINITY;
+        static constexpr bool has_quiet_NaN = FLT128_HAS_QUIET_NAN;
+        static constexpr bool has_signaling_NaN = has_quiet_NaN;
+        static constexpr float_denorm_style has_denorm = bool(FLT128_HAS_DENORM) ? denorm_present : denorm_absent;
+        static constexpr bool has_denorm_loss = std::numeric_limits<float>::has_denorm_loss;
+
+        static constexpr __float128 infinity()      { return std::numeric_limits<float>::infinity(); }
+        static constexpr __float128 quiet_NaN()     { return std::numeric_limits<float>::quiet_NaN(); }
+        static constexpr __float128 signaling_NaN() { return std::numeric_limits<float>::signaling_NaN(); }
+        static constexpr __float128 denorm_min()    { return FLT128_DENORM_MIN; }
+
+        static constexpr bool is_iec559 = has_infinity && has_quiet_NaN && has_denorm == denorm_present;
+        static constexpr bool is_bounded = true;
+        static constexpr bool is_modulo = false;
+
+        static constexpr bool traps = false;
+        static constexpr bool tinyness_before = false;
+        static constexpr float_round_style round_style = round_to_nearest;
+    };
+
+
+    inline __float128 abs(__float128 x) { return cabsq(__complex128{x,0.0Q}); }
+    inline __float128 acos(__float128 x) { return acosq(x); }
+    inline __float128 acosh(__float128 x) { return acoshq(x); }
+    inline __float128 asin(__float128 x) { return asinq(x); }
+    inline __float128 asinh(__float128 x) { return asinhq(x); }
+    inline __float128 atan(__float128 x) { return atanq(x); }
+    inline __float128 atanh(__float128 x) { return atanhq(x); }
+    inline __float128 cbrt(__float128 x) { return cbrtq(x); }
+    inline __float128 ceil(__float128 x) { return ceilq(x); }
+    inline __float128 cosh(__float128 x) { return coshq(x); }
+    inline __float128 cos(__float128 x) { return cosq(x); }
+    inline __float128 erf(__float128 x) { return erfq(x); }
+    inline __float128 erfc(__float128 x) { return erfcq(x); }
+    inline __float128 exp(__float128 x) { return expq(x); }
+    inline __float128 expm1(__float128 x) { return expm1q(x); }
+    inline __float128 fabs(__float128 x) { return fabsq(x); }
+    inline int        finite(__float128 x) { return finiteq(x); }
+    inline __float128 floor(__float128 x) { return floorq(x); }
+    inline __float128 frexp(__float128 x, int* p) { return frexpq(x,p); }
+    inline int        isinf(__float128 x) { return isinfq(x); }
+    inline int        ilogb(__float128 x) { return ilogbq(x); }
+    inline int        isnan(__float128 x) { return isnanq(x); }
+    inline __float128 j0(__float128 x) { return j0q(x); }
+    inline __float128 j1(__float128 x) { return j1q(x); }
+    inline __float128 jn(int i, __float128 x) { return jnq(i,x); }
+    inline __float128 ldexp(__float128 x, int i) { return ldexpq(x,i); }
+    inline __float128 lgamma(__float128 x) { return lgammaq(x); }
+    inline long long int llrint(__float128 x) { return llrintq(x); }
+    inline long long int llround(__float128 x) { return llroundq(x); }
+    inline __float128 log(__float128 x) { return logq(x); }
+    inline __float128 log10(__float128 x) { return log10q(x); }
+    inline __float128 log2(__float128 x) { return log2q(x); }
+    inline __float128 log1p(__float128 x) { return log1pq(x); }
+    inline long int   lrint(__float128 x) { return lrintq(x); }
+    inline long int   lround(__float128 x) { return lroundq(x); }
+    inline __float128 nearbyint(__float128 x) { return nearbyintq(x); }
+    inline __float128 pow(__float128 x, __float128 y) { return powq(x,y); }
+    inline __float128 rint(__float128 x) { return rintq(x); }
+    inline __float128 round(__float128 x) { return roundq(x); }
+    inline __float128 scalbln(__float128 x, long int li) { return scalblnq(x,li); }
+    inline __float128 scalbn(__float128 x, int i) { return scalbnq(x,i); }
+    inline int        signbit(__float128 x) { return signbitq(x); }
+    inline __float128 sinh(__float128 x) { return sinhq(x); }
+    inline __float128 sin(__float128 x) { return sinq(x); }
+    inline __float128 sqrt(__float128 x) { return sqrtq(x); }
+    inline __float128 tan(__float128 x) { return tanq(x); }
+    inline __float128 tanh(__float128 x) { return tanhq(x); }
+    inline __float128 tgamma(__float128 x) { return tgammaq(x); }
+    inline __float128 trunc(__float128 x) { return truncq(x); }
+    inline __float128 y0(__float128 x) { return y0q(x); }
+    inline __float128 y1(__float128 x) { return y1q(x); }
+    inline __float128 yn(int i, __float128 x) { return ynq(i,x); }
+
+
+    /* Prototypes for complex functions */
+    inline __float128 abs(__complex128 x) { return cabsq(x); }
+    inline __float128 arg(__complex128 x) { return cargq(x); }
+    inline __float128 imag(__complex128 x) { return cimagq(x); }
+    inline __float128 real(__complex128 x) { return crealq(x); }
+    inline __complex128 acos(__complex128 x) { return cacosq(x); }
+    inline __complex128 acosh(__complex128 x) { return cacoshq(x); }
+    inline __complex128 asin(__complex128 x) { return casinq(x); }
+    inline __complex128 asinh(__complex128 x) { return casinhq(x); }
+    inline __complex128 atan(__complex128 x) { return catanq(x); }
+    inline __complex128 atanh(__complex128 x) { return catanhq(x); }
+    inline __complex128 cos(__complex128 x) { return ccosq(x); }
+    inline __complex128 cosh(__complex128 x) { return ccoshq(x); }
+    inline __complex128 exp(__complex128 x) { return cexpq(x); }
+    inline __complex128 expi(__float128 x) { return cexpiq(x); }
+    inline __complex128 log10(__complex128 x) { return clog10q(x); }
+    inline __complex128 conj(__complex128 x) { return conjq(x); }
+    inline __complex128 pow(__complex128 x, __complex128 y) { return cpowq(x,y); }
+    inline __complex128 proj(__complex128 x) { return cprojq(x); }
+    inline __complex128 sin(__complex128 x) { return csinq(x); }
+    inline __complex128 sinh(__complex128 x) { return csinhq(x); }
+    inline __complex128 sqrt(__complex128 x) { return csqrtq(x); }
+    inline __complex128 tan(__complex128 x) { return ctanq(x); }
+    inline __complex128 tanh(__complex128 x) { return ctanhq(x); }
+
+    inline std::ostream& operator<<(std::ostream& os, __float128 x) {
+        const int prec = (os.precision()==0 ? 42 : os.precision());
+
+        char buf[128];
+        int n = quadmath_snprintf(buf, 128, "%+-#46.*Qe", prec, x);
+
+        if (n>128) {
+            char *str = new char[n+1];
+            if (str)
+                quadmath_snprintf (str, n+1, "%+-#46.*Qe", prec, x);
+            os << str;
+            delete(str);
+        }
+        else {
+            os << buf;
+        }
+        return os;
+    }
+}
+
+#include <complex>
+
+namespace std {
+
+    inline __float128 abs(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        return cabsq(X);
+    }
+    inline __float128 carg(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        return cargq(X);
+    }
+    inline __float128 cimag(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        return cimagq(X);
+    }
+    inline __float128 creal(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        return crealq(X);
+    }
+    inline std::complex<__float128> acos(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = cacosq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> acosh(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = cacoshq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> asin(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = casinq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> asinh(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = casinhq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> atan(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = catanq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> atanh(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = catanhq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> cos(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = ccosq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> cosh(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = ccoshq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> exp(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = cexpq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    //inline std::complex<__float128> expi(__float128 x) { 
+        //__complex128 X = cexpiq(x);
+        //return std::complex<__float128>(__real__ X, __imag__ X);
+    //}
+    inline std::complex<__float128> log10(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = clog10q(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> conj(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = conjq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> pow(std::complex<__float128> x, std::complex<__float128> y) {
+        __complex128 X{x.real(),x.imag()};
+        __complex128 Y{y.real(),y.imag()};
+        X = cpowq(X,Y);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> proj(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = cprojq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> sin(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = csinq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> sinh(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = csinhq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> sqrt(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = csqrtq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> tan(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = ctanq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+    inline std::complex<__float128> tanh(std::complex<__float128> x) { 
+        __complex128 X{x.real(),x.imag()};
+        X = ctanhq(X);
+        return std::complex<__float128>(__real__ X, __imag__ X);
+    }
+
+    inline std::complex< __float128 > pow(std::complex< __float128> x , __float128 y) {
+        __float128 R = powq(abs(x), y);
+        __float128 phi = atanq(x.imag()/x.real());
+        return std::complex<__float128 >(R*cosq(y*phi), R*sinq(y*phi));
+    }
+}
+
+#endif /* end of include guard: HYSOP_QUAD_MATHS_H */
+
+#endif
diff --git a/HySoP/src/hysop++/src/utils/utils.h b/HySoP/src/hysop++/src/utils/utils.h
index 066516c2836e611aa2db93072dc3601a6feb45a4..ead393fdfcbeec120608a0d9da2a161a0a9fddb2 100644
--- a/HySoP/src/hysop++/src/utils/utils.h
+++ b/HySoP/src/hysop++/src/utils/utils.h
@@ -9,6 +9,7 @@
 #include <iostream>
 #include <limits>
 
+#include "maths/quad_maths.h"
 #include "detail/index_seq.h"
 #include "boost/multi_array.hpp"