Merge pull request #10 from itpplasma/7-make-mfem-optional-dependency

Make MFEM optional dependency

Author: krystophny

CMakeLists.txt

@@ -5,6 +5,9 @@ set(CMAKE_MODULE_PATH ${CMAKE_MODULE_PATH} ${CMAKE_SOURCE_DIR}/cmake)
 ### Define the project
 project(MEPHIT LANGUAGES C CXX Fortran)
 
+### Define options
+option(WITH_MFEM "Compile with MFEM" FALSE)
+
 ### Specify paths
 set(EXECUTABLE_OUTPUT_PATH ${PROJECT_BINARY_DIR}/bin)
 set(LIBRARY_OUTPUT_PATH ${PROJECT_BINARY_DIR}/lib)
@@ -90,7 +93,10 @@ add_link_options(${NETCDF_FLIBS})
 include(ExternalProject)
 include(SetupCODE)
 include(SetupTriangle)
-include(SetupMFEM)
+if(WITH_MFEM)
+  include(SetupMFEM)
+  add_definitions(-DUSE_MFEM)
+endif()
 include(Util)
 find_or_fetch(libneo)
 
@@ -128,9 +134,9 @@ set_source_files_properties(${COMMON_FORTRAN_SRC_FILES}
 set(MEPHIT_FORTRAN_SRC_FILES
   src/mephit_conf.f90
   src/mephit_util.f90
-  src/mephit_mesh.f90
+  src/mephit_mesh.F90
   src/mephit_pert.f90
-  src/mephit_iter.f90
+  src/mephit_iter.F90
 )
 set_source_files_properties(src/mephit_test.f90 ${MEPHIT_FORTRAN_SRC_FILES}
   PROPERTIES COMPILE_FLAGS "${CMAKE_Fortran_FLAGS} ${WARNING_FLAGS}")
@@ -142,9 +148,13 @@ set(MEPHIT_C_SRC_FILES
 )
 set(MEPHIT_CPP_SRC_FILES
   src/mephit_fem.cpp
-  src/magnetic_differential_equation.cpp
-  src/miscellaneous.cpp
 )
+if(WITH_MFEM)
+  set(MEPHIT_CPP_SRC_FILES ${MEPHIT_CPP_SRC_FILES}
+    src/magnetic_differential_equation.cpp
+    src/miscellaneous.cpp
+  )
+endif()
 set_source_files_properties(src/mephit_run.c ${MEPHIT_C_SRC_FILES} ${MEPHIT_CPP_SRC_FILES}
 	PROPERTIES COMPILE_FLAGS "${CMAKE_C_FLAGS} -DREAL=double -I${TRIANGLE_INCLUDE_DIR} -I${SUITESPARSE_INCLUDE_DIRS} -L${TRIANGLE_LIB_DIR}")
 
@@ -157,11 +167,10 @@ add_library(mephit SHARED
   ${MEPHIT_C_SRC_FILES}
   ${MEPHIT_CPP_SRC_FILES}
 )
-add_dependencies(mephit MFEM TRIANGLE)
+add_dependencies(mephit TRIANGLE)
 target_include_directories(mephit PUBLIC
   ${magfie_include_dir}
   ${NETCDF_INCLUDES}
-  ${MFEM_INCLUDE_DIRS}
   ${NETCDF_INCLUDES}
   ${TRIANGLE_INCLUDE_DIR}
 )
@@ -171,7 +180,6 @@ set(MEPHIT_LIBS
   hdf5::hdf5_fortran
   hdf5::hdf5_hl_fortran
   GSL::gsl
-  ${MFEM_LIBRARIES}
   BLAS::BLAS
   LAPACK::LAPACK
   ${FFTW_LIBRARIES}
@@ -186,6 +194,11 @@ if(LEGACY_SUITESPARSE)
 else()
   target_link_libraries(mephit SuiteSparse::umfpack)
 endif()
+if(WITH_MFEM)
+  add_dependencies(mephit MFEM)
+  target_include_directories(mephit PUBLIC ${MFEM_INCLUDE_DIRS})
+  target_link_libraries(mephit ${MFEM_LIBRARIES})
+endif()
 
 ### Define executables
 add_executable(mephit_run.x src/mephit_run.c)

README.md

@@ -15,7 +15,7 @@ Prerequisites from external sources for running MEPHIT are as follows.
 - [Triangle](https://www.cs.cmu.edu/~quake/triangle.html)
 - [Boost](https://www.boost.org/)
 - [FreeFem++](https://github.com/FreeFem/FreeFem-sources)
-- [MFEM](https://mfem.org/)
+- [MFEM](https://mfem.org/) is optional
 - [HDF5](https://www.hdfgroup.org/downloads/hdf5)
 - [NetCDF](https://github.com/Unidata/netcdf-fortran)
 - Python 3 including packages listed in [`requirements.txt`](requirements.txt) for plotting

src/mephit_conf.f90

@@ -198,9 +198,6 @@ module mephit_conf
     !> Enable debugging of initial iterations without plasma response. Defaults to true.
     logical :: debug_initial = .true.
 
-    !> Enable MFEM interface. Defaults to false.
-    logical :: debug_mfem = .false.
-
     !> Single poloidal mode used in comparison with KiLCA code. Defaults to 0 (ASDEX
     !> geometry).
     integer :: kilca_pol_mode = 0
@@ -348,8 +345,6 @@ subroutine config_export_hdf5(config, file, group)
       comment = 'offset added to radial electric field', unit = 'statV cm^-1')
     call h5_add(h5id_root, grp // '/debug_initial', config%debug_initial, &
       comment = 'enable debugging of initial runs without plasma response')
-    call h5_add(h5id_root, grp // '/debug_MFEM', config%debug_MFEM, &
-      comment = 'enable MFEM interface')
     call h5_add(h5id_root, grp // '/kilca_pol_mode', config%kilca_pol_mode, &
       comment = 'single poloidal mode used in comparison with KiLCA code')
     call h5_add(h5id_root, grp // '/kilca_scale_factor', config%kilca_scale_factor, &

src/mephit_fem.cpp

@@ -1,6 +1,8 @@
 #include "mephit_fem.h"
+#ifdef USE_MFEM
 #include "mfem.hpp"
 #include "magnetic_differential_equation.h"
+#endif  // USE_MFEM
 #include <boost/geometry.hpp>
 #include <boost/geometry/geometries/point_xy.hpp>
 #include <boost/geometry/index/rtree.hpp>
@@ -39,6 +41,8 @@ extern "C" void Rtree_query(double R, double Z, int *result_size, int **results)
   *results = query_results.data();
 }
 
+#ifdef USE_MFEM
+
 typedef std::map<std::pair<double, double>, size_t> points_2D;
 
 /*
@@ -142,3 +146,5 @@ extern "C" int FEM_test(const char *mesh_file,
   }
   return 0;
 }
+
+#endif  // USE_MFEM

src/mephit_fem.h

@@ -50,12 +50,15 @@ typedef void real_vector_field(const double R,
 typedef void complex_scalar_field(const double R,
                                   const double Z,
                                   complex_double *scalar);
+
+#ifdef USE_MFEM
 int FEM_test(const char *mesh_file,
              const int tor_mode,
              const int n_dof,
              complex_double *dof,
              real_vector_field *unit_B0,
              complex_scalar_field *MDE_inhom);
+#endif  // USE_MFEM
 
 #ifdef __cplusplus
 }

src/mephit_iter.f90 renamed to src/mephit_iter.F90

@@ -541,12 +541,14 @@ subroutine perteq_iterate(perteq, precond, fdm)
       write (postfix, postfix_fmt) kiter
       Bn_prev%DOF(:) = perteq%Bn%DOF
       Bn_prev%comp_phi(:) = perteq%Bn%comp_phi
-      ! compute B_(n+1) = K * B_n + B_vac ... different from next_iteration_arnoldi
-      if (kiter <= 1 .and. conf%debug_MFEM) then
+#ifdef USE_MFEM
+      if (kiter <= 1) then
         call MFEM_test(perteq%pn)
         call perteq_write('("iter/", a, "MFEM_' // postfix // '")', &
             ' (after MFEM iteration)', presn = perteq%pn, presmn = perteq%pn)
       end if
+#endif
+      ! compute B_(n+1) = K * B_n + B_vac ... different from next_iteration_arnoldi
       call compute_presn(perteq, fdm, conf%damp)
       if (kiter <= 1) then
         call perteq_write('("iter/", a, "' // postfix // '")', &
@@ -634,16 +636,16 @@ subroutine debug_initial_iteration(perteq, fdm)
     if (conf%debug_initial) then
       perteq%Bn%DOF(:) = vac%Bn%DOF
       perteq%Bn%comp_phi(:) = vac%Bn%comp_phi
-      if (conf%debug_MFEM) then
-        call MFEM_test(perteq%pn)
-        call perteq_write('("debug_MFEM_initial/MFEM_", a)', &
-          ' (initial MFEM iteration)', presn = perteq%pn, presmn = perteq%pn)
-      end if
+#ifdef USE_MFEM
+      call MFEM_test(perteq%pn)
+      call perteq_write('("debug_MFEM_initial/MFEM_", a)', &
+        ' (initial MFEM iteration)', presn = perteq%pn, presmn = perteq%pn)
+#endif
       call compute_presn(perteq, fdm, .false.)
-      if (conf%debug_MFEM) then
-        call perteq_write('("debug_MFEM_initial/", a)', &
-          ' (initial iteration)', presn = perteq%pn, presmn = perteq%pn)
-      end if
+#ifdef USE_MFEM
+      call perteq_write('("debug_MFEM_initial/", a)', &
+        ' (initial iteration)', presn = perteq%pn, presmn = perteq%pn)
+#endif
       call compute_currn(perteq, fdm, .false., .true.)
       perteq%Bn%DOF(:) = vac%Bn%DOF
       perteq%Bn%comp_phi(:) = vac%Bn%comp_phi
@@ -705,6 +707,7 @@ subroutine presn_inhom(R, Z, scalar) bind(C, name = 'presn_inhom')
     scalar = -dp0_dpsi * (B_n(1) * B_0(3) - B_n(3) * B_0(1)) * R / sqrt(sum(B_0 * B_0))
   end subroutine presn_inhom
 
+#ifdef USE_MFEM
   subroutine MFEM_test(pn)
     use iso_c_binding, only: c_int, c_null_char, c_loc, c_funloc
     use mephit_conf, only: conf, logger, basename_suffix, decorate_filename
@@ -720,6 +723,7 @@ subroutine MFEM_test(pn)
       call logger%write_msg
     end if
   end subroutine MFEM_test
+#endif
 
   subroutine FDM_init(fdm, nnz)
     type(FDM_t), intent(inout) :: fdm

src/mephit_mesh.f90 renamed to src/mephit_mesh.F90

@@ -1004,9 +1004,9 @@ subroutine generate_mesh
     call init_flux_variables
     call compare_gpec_coordinates
     call connect_mesh_points
-    if (conf%debug_MFEM) then
-      call mesh_write_MFEM
-    end if
+#ifdef USE_MFEM
+    call mesh_write_MFEM
+#endif
     call write_FreeFem_mesh
     call cache_init(cache, 4)
     call cache_equilibrium_field
@@ -2843,6 +2843,7 @@ subroutine mesh_write(mesh, file, dataset)
     call h5_close(h5id_root)
   end subroutine mesh_write
 
+#ifdef USE_MFEM
   subroutine mesh_write_MFEM
     use mephit_conf, only: basename_suffix, decorate_filename
     integer :: fid, ktri, kp, kpoi
@@ -2870,6 +2871,7 @@ subroutine mesh_write_MFEM
     end do
     close(fid)
   end subroutine mesh_write_MFEM
+#endif
 
   subroutine write_FreeFem_mesh
     use iso_c_binding, only: c_null_char