diff --git a/CMakeLists.txt b/CMakeLists.txt
index f6865bba8194db58c889e2e17e81c8cbbeaabc5c..047b0f97cba73a12d3091877931f771f2385a3c8 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -108,8 +108,6 @@ if (OPT_ENABLE_GPU_SOLVER)
## Find HIP stuff
find_package(HIP QUIET MODULE)
if (HIP_FOUND)
- execute_process(COMMAND sh -c "/opt/rocm/hip/bin/hipify-perl ${CMAKE_SOURCE_DIR}/src/kernels_cuda.cu > ${CMAKE_SOURCE_DIR}/src/hipified_kernels_cuda.cpp")
- execute_process(COMMAND sh -c "/opt/rocm/hip/bin/hipify-perl ${CMAKE_SOURCE_DIR}/src/maxwell_kernels_cuda.cu > ${CMAKE_SOURCE_DIR}/src/hipified_maxwell_kernels_cuda.cpp")
add_definitions(-DHAVE_HIP)
include_directories("/opt/rocm/include")
set(CMAKE_HIP_LINK_EXECUTABLE "${HIP_HIPCC_CMAKE_LINKER_HELPER} ${HCC_HOME} <FLAGS> <CMAKE_CXX_LINK_FLAGS> <LINK_FLAGS> <OBJECTS> -o <TARGET> <LINK_LIBRARIES>")
@@ -123,6 +121,8 @@ if (OPT_ENABLE_GPU_SOLVER)
set(HIP_NVCC_FLAGS "-ccbin g++-8" ${HIP_NVCC_FLAGS})
option(OPT_USE_HIP "Use HIP GPU code path" OFF)
if (OPT_USE_HIP)
+ execute_process(COMMAND sh -c "/opt/rocm/hip/bin/hipify-perl ${CMAKE_SOURCE_DIR}/src/kernels_cuda.cu > ${CMAKE_SOURCE_DIR}/src/hipified_kernels_cuda.cpp")
+ execute_process(COMMAND sh -c "/opt/rocm/hip/bin/hipify-perl ${CMAKE_SOURCE_DIR}/src/maxwell_kernels_cuda.cu > ${CMAKE_SOURCE_DIR}/src/hipified_maxwell_kernels_cuda.cpp")
add_definitions(-DGPU_USE_HIP)
add_definitions(-D__HIP_PLATFORM_HCC__)
set(LINK_LIBS ${LINK_LIBS})
@@ -163,8 +163,8 @@ if (COMPILER_IS_CLANG)
-Wno-padded -Wno-shorten-64-to-32 -Wno-sign-conversion -Wno-old-style-cast \
-Wno-sign-compare -Wno-c99-extensions -Wno-extra-semi-stmt -Wno-source-uses-openmp")
else()
- set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wall -Wno-sign-compare -Wshadow \
- -Wno-unknown-pragmas")
+ set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -Wall -Wextra -Wno-sign-compare -Wshadow \
+ -Wno-unknown-pragmas -Wno-unused-parameter")
endif()
set(CMAKE_CXX_FLAGS_DEBUG "-O0 -g -fpermissive")
set(CMAKE_CXX_FLAGS_RELEASE "-O3 -g -DNDEBUG -fpermissive")
diff --git a/doc/changelog.md b/doc/changelog.md
index cc1ee8067da36fa2b702fa49a05435a0330f1f4f..ccdb90e9d1f0cbe0a9a06bc3989ee25770924e84 100644
--- a/doc/changelog.md
+++ b/doc/changelog.md
@@ -18,5 +18,5 @@ Unavailable features:
- Switched back to C++17 due to various problems on Eigen and Mac OS X.
- Use system GMSH if installed and only after look at the hardcoded path (a55af8ab).
- Improved export to SILO. E, H, and J fields are now exportable as nodal or zonal variables. It is also possible to disable the export of a specific field, see `lua_api.md` (8a831181).
-- Preliminary implementation of leapfrog integration on CPU
+- Preliminary implementation of leapfrog integration on CPU.
diff --git a/doc/lua_api.md b/doc/lua_api.md
index 8610533dd2b8308a62ed5bd6d8e5149e57cbd1a3..363323459c11af48f402e5a8a34edebdf34392b3 100644
--- a/doc/lua_api.md
+++ b/doc/lua_api.md
@@ -20,7 +20,7 @@ This document describes the API available on the version v0.1 of the solver.
- `sim.use_gpu` (0/1): enable/disable GPU usage.
- `sim.approx_order` (integer): method approximation order.
- `sim.geom_order` (integer): geometry order. Only order 1 supported for now.
-- sim.time_integrator (string): `"rk4"`, `"leapfrog"` or `"euler"`. Euler is only for test purporses.
+- `sim.time_integrator` (string): `"rk4"`, `"leapfrog"` or `"euler"`. Euler is only for test purporses, as it is numerically instable. Leapfrog works on CPU for now, use it if you have dissipative materials.
### Postprocessing general variables
diff --git a/include/param_loader.h b/include/param_loader.h
index f656a71a596586e162e481bd1c74f9555b9cd4a6..fd2f5eac9b7af2da537533d2f98c0286e57dc16c 100644
--- a/include/param_loader.h
+++ b/include/param_loader.h
@@ -46,6 +46,7 @@ public:
std::string sim_gmshmodel(void) const;
bool sim_usegpu(void) const;
time_integrator_type sim_timeIntegrator(void) const;
+ std::string sim_timeIntegratorName(void) const;
size_t postpro_siloOutputRate(void) const;
size_t postpro_cyclePrintRate(void) const;
diff --git a/src/maxwell_cpu.cpp b/src/maxwell_cpu.cpp
index 35d8e642b569f9299d6688aefe6f26c46e2506e0..a5a13a54949a64bf9298d694bc045d03b5a9cb48 100644
--- a/src/maxwell_cpu.cpp
+++ b/src/maxwell_cpu.cpp
@@ -287,7 +287,6 @@ compute_fluxes_planar(solver_state& state)
}
}
-/* In leapfrog we use centered fluxes for the moment */
static void
compute_fluxes_planar_E(solver_state& state)
{
@@ -310,22 +309,27 @@ compute_fluxes_planar_E(solver_state& state)
auto gofs = ed.flux_base + lofs;
/* Sources are imposed on jumps */
+ auto jEx = state.jumps.Ex[gofs] - state.bndsrcs.Ex[gofs];
+ auto jEy = state.jumps.Ey[gofs] - state.bndsrcs.Ey[gofs];
+ auto jEz = state.jumps.Ez[gofs] - state.bndsrcs.Ez[gofs];
auto jHx = state.jumps.Hx[gofs] + state.bndsrcs.Hx[gofs];
auto jHy = state.jumps.Hy[gofs] + state.bndsrcs.Hy[gofs];
auto jHz = state.jumps.Hz[gofs] + state.bndsrcs.Hz[gofs];
- auto aE = face_det * 0.5;
+ auto ndotE = nx*jEx + ny*jEy + nz*jEz;
+
+ auto aE = face_det * state.matparams.aE[gofs];
+ auto bE = face_det * state.matparams.bE[gofs];
/* Compute fluxes */
- state.fluxes.Ex[gofs] = aE*(nz*jHy - ny*jHz);
- state.fluxes.Ey[gofs] = aE*(nx*jHz - nz*jHx);
- state.fluxes.Ez[gofs] = aE*(ny*jHx - nx*jHy);
+ state.fluxes.Ex[gofs] = aE*(nz*jHy - ny*jHz) + bE*(ndotE*nx - jEx);
+ state.fluxes.Ey[gofs] = aE*(nx*jHz - nz*jHx) + bE*(ndotE*ny - jEy);
+ state.fluxes.Ez[gofs] = aE*(ny*jHx - nx*jHy) + bE*(ndotE*nz - jEz);
}
}
}
}
-/* In leapfrog we use centered fluxes for the moment */
static void
compute_fluxes_planar_H(solver_state& state)
{
@@ -351,13 +355,19 @@ compute_fluxes_planar_H(solver_state& state)
auto jEx = state.jumps.Ex[gofs] - state.bndsrcs.Ex[gofs];
auto jEy = state.jumps.Ey[gofs] - state.bndsrcs.Ey[gofs];
auto jEz = state.jumps.Ez[gofs] - state.bndsrcs.Ez[gofs];
+ auto jHx = state.jumps.Hx[gofs] + state.bndsrcs.Hx[gofs];
+ auto jHy = state.jumps.Hy[gofs] + state.bndsrcs.Hy[gofs];
+ auto jHz = state.jumps.Hz[gofs] + state.bndsrcs.Hz[gofs];
+
+ auto ndotH = nx*jHx + ny*jHy + nz*jHz;
- auto aH = face_det * 0.5;
+ auto aH = face_det * state.matparams.aH[gofs];
+ auto bH = face_det * state.matparams.bH[gofs];
/* Compute fluxes */
- state.fluxes.Hx[gofs] = aH*(ny*jEz - nz*jEy);
- state.fluxes.Hy[gofs] = aH*(nz*jEx - nx*jEz);
- state.fluxes.Hz[gofs] = aH*(nx*jEy - ny*jEx);
+ state.fluxes.Hx[gofs] = aH*(ny*jEz - nz*jEy) + bH*(ndotH*nx - jHx);
+ state.fluxes.Hy[gofs] = aH*(nz*jEx - nx*jEz) + bH*(ndotH*ny - jHy);
+ state.fluxes.Hz[gofs] = aH*(nx*jEy - ny*jEx) + bH*(ndotH*nz - jHz);
}
}
}
diff --git a/src/maxwell_solver.cpp b/src/maxwell_solver.cpp
index 5ab9afd6b711ed93e1633c44b13c8d4afe6cb85d..69ee2b95ce30a2e91c9bebd66be8b562e1874b52 100644
--- a/src/maxwell_solver.cpp
+++ b/src/maxwell_solver.cpp
@@ -75,12 +75,16 @@ void test_it(const model& mod, State& state, maxwell::parameter_loader& mpl)
auto num_timesteps = mpl.sim_timesteps();
auto silo_output_rate = mpl.postpro_siloOutputRate();
auto cycle_print_rate = mpl.postpro_cyclePrintRate();
+ auto ti = mpl.sim_timeIntegrator();
#ifdef _OPENMP
omp_set_num_threads(4);
#endif
- auto ti = mpl.sim_timeIntegrator();
+ std::cout << " BEGINNING SIMULATION" << std::endl;
+ std::cout << "I will do " << num_timesteps << " of " << state.delta_t;
+ std::cout << "s each." << std::endl;
+ std::cout << "Time integrator: " << mpl.sim_timeIntegratorName() << std::endl;
prepare_sources(mod, state, mpl);
for(size_t i = 0; i < num_timesteps; i++)
diff --git a/src/param_loader.cpp b/src/param_loader.cpp
index 0af01bd73006d123e9b02ee952510872dd6b6026..b89bf4f963f5b30f57019bfbc527e99dfd30a070 100644
--- a/src/param_loader.cpp
+++ b/src/param_loader.cpp
@@ -233,4 +233,12 @@ parameter_loader_base::sim_timeIntegrator(void) const
std::cout << "[CONFIG] warning: invalid time integrator '";
std::cout << ti << "'" << std::endl;
return time_integrator_type::RK4;
-}
\ No newline at end of file
+}
+
+std::string
+parameter_loader_base::sim_timeIntegratorName(void) const
+{
+ std::string ti = lua["sim"]["time_integrator"];
+ return ti;
+}
+