diff --git a/core/docs/Input.md b/core/docs/Input.md index 6a5c8d829..a8693ecb2 100644 --- a/core/docs/Input.md +++ b/core/docs/Input.md @@ -185,6 +185,8 @@ ddi_n_periodic_images 4 4 4 ### DDI cutoff radius (if cutoff is used) ddi_radius 0.0 + +ddi_pb_zero_padding 1.0 ``` *Anisotropy:* @@ -197,7 +199,7 @@ a magnitude `K`. Via the keyword `ddi_method` the method employed to calculate the dipole-dipole interactions is specified. `none` - Dipole-Dipole interactions are neglected - `fft` - Uses a fast convolution method to accelerate the calculation + `fft` - Uses a fast convolution method to accelerate the calculation (RECOMMENDED) `cutoff` - Lets only spins within a maximal distance of 'ddi_radius' interact `fmm` - Uses the Fast-Multipole-Method (NOT YET IMPLEMENTED!) @@ -208,6 +210,10 @@ If the boundary conditions are periodic `ddi_n_periodic_images` specifies how ma *Note:* The images are appended on both sides (the edges get filled too) i.e. 1 0 0 -> one image in +a direction and one image in -a direction +If the boundary conditions are open in a lattice direction and sufficiently many periodic images are chosen, zero-padding in that direction can be skipped. +This improves the speed and memory footprint of the calculation, but comes at the cost of a very slight asymmetry in the interactions (decreasing with increasing periodic images). +If `ddi_pb_zero_padding` is set to 1, zero-padding is performed - even if the boundary condition is periodic in a direction. If it is set to 0, zero-padding is skipped. + **Neighbour shells**: Using `hamiltonian heisenberg_neighbours`, pair-wise interactions are handled in terms of diff --git a/core/include/Spirit/Hamiltonian.h b/core/include/Spirit/Hamiltonian.h index 450158609..98bf8729b 100644 --- a/core/include/Spirit/Hamiltonian.h +++ b/core/include/Spirit/Hamiltonian.h @@ -85,8 +85,9 @@ Configure the dipole-dipole interaction boundary conditions are used - `cutoff_radius`: the distance at which to stop the direct summation, if used +- `pb_zero_padding`: if `True` zero padding is used even for periodical directions */ -PREFIX void Hamiltonian_Set_DDI(State *state, int ddi_method, int n_periodic_images[3], float cutoff_radius=0, int idx_image=-1, int idx_chain=-1) SUFFIX; +PREFIX void Hamiltonian_Set_DDI(State *state, int ddi_method, int n_periodic_images[3], float cutoff_radius=0, bool pb_zero_padding=true, int idx_image=-1, int idx_chain=-1) SUFFIX; /* Getters @@ -133,13 +134,12 @@ PREFIX int Hamiltonian_Get_DMI_N_Pairs(State *state, int idx_image=-1, int idx_ Retrieves the dipole-dipole interaction configuration. - `ddi_method`: see integers defined above -- `n_periodic_images`: how many repetition of the spin configuration to - append along the translation directions [a, b, c], if periodical - boundary conditions are used -- `cutoff_radius`: the distance at which to stop the direct summation, - if used +- `n_periodic_images`: how many repetitions of the spin configuration to + append along the translation directions [a, b, c], if periodical boundary conditions are used +- `cutoff_radius`: the distance at which to stop the direct summation, if method_cutoff is used +- `pb_zero_padding`: if `True` zero padding is used even for periodical directions */ -PREFIX void Hamiltonian_Get_DDI(State *state, int * ddi_method, int n_periodic_images[3], float * cutoff_radius, int idx_image=-1, int idx_chain=-1) SUFFIX; +PREFIX void Hamiltonian_Get_DDI(State *state, int * ddi_method, int n_periodic_images[3], float * cutoff_radius, bool * pb_zero_padding, int idx_image=-1, int idx_chain=-1) SUFFIX; #include "DLL_Undefine_Export.h" #endif \ No newline at end of file diff --git a/core/include/engine/Hamiltonian_Heisenberg.hpp b/core/include/engine/Hamiltonian_Heisenberg.hpp index b508156c0..0e9b14c27 100644 --- a/core/include/engine/Hamiltonian_Heisenberg.hpp +++ b/core/include/engine/Hamiltonian_Heisenberg.hpp @@ -35,7 +35,7 @@ namespace Engine intfield anisotropy_indices, scalarfield anisotropy_magnitudes, vectorfield anisotropy_normals, pairfield exchange_pairs, scalarfield exchange_magnitudes, pairfield dmi_pairs, scalarfield dmi_magnitudes, vectorfield dmi_normals, - DDI_Method ddi_method, intfield ddi_n_periodic_images, scalar ddi_radius, + DDI_Method ddi_method, intfield ddi_n_periodic_images, bool ddi_pb_zero_padding, scalar ddi_radius, quadrupletfield quadruplets, scalarfield quadruplet_magnitudes, std::shared_ptr geometry, intfield boundary_conditions @@ -46,7 +46,7 @@ namespace Engine intfield anisotropy_indices, scalarfield anisotropy_magnitudes, vectorfield anisotropy_normals, scalarfield exchange_shell_magnitudes, scalarfield dmi_shell_magnitudes, int dm_chirality, - DDI_Method ddi_method, intfield ddi_n_periodic_images, scalar ddi_radius, + DDI_Method ddi_method, intfield ddi_n_periodic_images, bool ddi_pb_zero_padding, scalar ddi_radius, quadrupletfield quadruplets, scalarfield quadruplet_magnitudes, std::shared_ptr geometry, intfield boundary_conditions @@ -100,6 +100,7 @@ namespace Engine // Dipole Dipole interaction DDI_Method ddi_method; intfield ddi_n_periodic_images; + bool ddi_pb_zero_padding; // ddi cutoff variables scalar ddi_cutoff_radius; pairfield ddi_pairs; @@ -177,12 +178,12 @@ namespace Engine FFT::StrideContainer spin_stride; FFT::StrideContainer dipole_stride; - //Calculate the FT of the padded D matriess + // Calculate the FT of the padded D-matrics void FFT_Dipole_Matrices(FFT::FFT_Plan & fft_plan_dipole, int img_a, int img_b, int img_c); - //Calculate the FT of the padded spins + // Calculate the FT of the padded spins void FFT_Spins(const vectorfield & spins); - //Bounds for nested for loops. Only important for the CUDA version + // Bounds for nested for loops. Only important for the CUDA version field it_bounds_pointwise_mult; field it_bounds_write_gradients; field it_bounds_write_spins; diff --git a/core/python/spirit/hamiltonian.py b/core/python/spirit/hamiltonian.py index 4c1c0fe3d..ebc5210c1 100644 --- a/core/python/spirit/hamiltonian.py +++ b/core/python/spirit/hamiltonian.py @@ -93,19 +93,20 @@ def set_dmi(p_state, n_shells, D_ij, chirality=CHIRALITY_BLOCH, idx_image=-1, id ctypes.c_int(chirality), ctypes.c_int(idx_image), ctypes.c_int(idx_chain)) _Set_DDI = _spirit.Hamiltonian_Set_DDI -_Set_DDI.argtypes = [ctypes.c_void_p, ctypes.c_int, ctypes.POINTER(ctypes.c_int), ctypes.c_float, +_Set_DDI.argtypes = [ctypes.c_void_p, ctypes.c_int, ctypes.POINTER(ctypes.c_int), ctypes.c_float, ctypes.c_bool, ctypes.c_int, ctypes.c_int] _Set_DDI.restype = None -def set_ddi(p_state, ddi_method, n_periodic_images=[4,4,4], radius=0.0, idx_image=-1, idx_chain=-1): +def set_ddi(p_state, ddi_method, n_periodic_images=[4,4,4], radius=0.0, pb_zero_padding=True, idx_image=-1, idx_chain=-1): """Set the dipolar interaction calculation method. - `ddi_method`: one of the integers defined above - `n_periodic_images`: the number of periodical images in the three translation directions, taken into account when boundaries in the corresponding direction are periodical - `radius`: the cutoff radius for the direct summation method + - `pb_zero_padding`: if `True` zero padding is used for periodical directions """ vec3 = ctypes.c_int * 3 - _Set_DDI(ctypes.c_void_p(p_state), ctypes.c_int(ddi_method) , vec3(*n_periodic_images), ctypes.c_float(radius), + _Set_DDI(ctypes.c_void_p(p_state), ctypes.c_int(ddi_method) , vec3(*n_periodic_images), ctypes.c_float(radius), ctypes.c_bool(pb_zero_padding), ctypes.c_int(idx_image), ctypes.c_int(idx_chain)) ### ---------------------------------- Get ---------------------------------- @@ -145,9 +146,16 @@ def get_field(p_state, idx_image=-1, idx_chain=-1): return float(magnitude), [n for n in normal] _Get_DDI = _spirit.Hamiltonian_Get_DDI -_Get_DDI.argtypes = [ctypes.c_void_p, ctypes.c_int, ctypes.c_int] -_Get_DDI.restype = ctypes.c_float +_Get_DDI.argtypes = [ctypes.c_void_p, ctypes.POINTER(ctypes.c_int), ctypes.POINTER(ctypes.c_int), ctypes.POINTER(ctypes.c_float), ctypes.POINTER(ctypes.c_bool), ctypes.c_int, ctypes.c_int] +_Get_DDI.restype = None def get_ddi(p_state, idx_image=-1, idx_chain=-1): - """Returns the cutoff radius of the DDI.""" - return float(_Get_DDI(ctypes.c_void_p(p_state), ctypes.c_int(idx_image), - ctypes.c_int(idx_chain))) \ No newline at end of file + """Returns a dictionary, containing information about the used ddi settings""" + ddi_method = ctypes.c_int() + n_periodic_images = (3*ctypes.c_int)() + cutoff_radius = ctypes.c_float() + pb_zero_padding = ctypes.c_bool() + _Get_DDI(ctypes.c_void_p(p_state), ctypes.byref(ddi_method), n_periodic_images, ctypes.byref(cutoff_radius), ctypes.byref(pb_zero_padding), ctypes.c_int(idx_image), ctypes.c_int(idx_chain)) + return { "ddi_method" : ddi_method.value, + "n_periodic_images" : [ i for i in n_periodic_images ], + "cutoff_radius" : cutoff_radius.value, + "pb_zero_padding" : pb_zero_padding.value } \ No newline at end of file diff --git a/core/src/Spirit/Hamiltonian.cpp b/core/src/Spirit/Hamiltonian.cpp index 9955e6b14..446c2ee1b 100644 --- a/core/src/Spirit/Hamiltonian.cpp +++ b/core/src/Spirit/Hamiltonian.cpp @@ -267,7 +267,7 @@ void Hamiltonian_Set_DMI(State *state, int n_shells, const float * dij, int chir } } -void Hamiltonian_Set_DDI(State *state, int ddi_method, int n_periodic_images[3], float cutoff_radius, int idx_image, int idx_chain) noexcept +void Hamiltonian_Set_DDI(State *state, int ddi_method, int n_periodic_images[3], float cutoff_radius, bool pb_zero_padding, int idx_image, int idx_chain) noexcept { try { @@ -289,11 +289,12 @@ void Hamiltonian_Set_DDI(State *state, int ddi_method, int n_periodic_images[3], ham->ddi_n_periodic_images[1] = n_periodic_images[1]; ham->ddi_n_periodic_images[2] = n_periodic_images[2]; ham->ddi_cutoff_radius = cutoff_radius; + ham->ddi_pb_zero_padding = pb_zero_padding; ham->Update_Interactions(); Log( Utility::Log_Level::Info, Utility::Log_Sender::API, fmt::format( - "Set ddi to method {}, periodic images {} {} {} and cutoff radius {}", - ddi_method, n_periodic_images[0], n_periodic_images[1], n_periodic_images[2], cutoff_radius), idx_image, idx_chain ); + "Set ddi to method {}, periodic images {} {} {}, cutoff radius {} and pb_zero_padding {}", + ddi_method, n_periodic_images[0], n_periodic_images[1], n_periodic_images[2], cutoff_radius, pb_zero_padding), idx_image, idx_chain ); } else Log( Utility::Log_Level::Warning, Utility::Log_Sender::API, "DDI cannot be set on " + @@ -552,7 +553,7 @@ int Hamiltonian_Get_DMI_N_Pairs(State *state, int idx_image, int idx_chain) noe } } -void Hamiltonian_Get_DDI(State *state, int * ddi_method, int n_periodic_images[3], float * cutoff_radius, int idx_image, int idx_chain) noexcept +void Hamiltonian_Get_DDI(State *state, int * ddi_method, int n_periodic_images[3], float * cutoff_radius, bool * pb_zero_padding, int idx_image, int idx_chain) noexcept { try { @@ -571,6 +572,7 @@ void Hamiltonian_Get_DDI(State *state, int * ddi_method, int n_periodic_images[3 n_periodic_images[1] = (int)ham->ddi_n_periodic_images[1]; n_periodic_images[2] = (int)ham->ddi_n_periodic_images[2]; *cutoff_radius = (float)ham->ddi_cutoff_radius; + *pb_zero_padding = ham->ddi_pb_zero_padding; } } catch( ... ) diff --git a/core/src/engine/Hamiltonian_Heisenberg.cpp b/core/src/engine/Hamiltonian_Heisenberg.cpp index d6c1f1dea..be1d9a71b 100644 --- a/core/src/engine/Hamiltonian_Heisenberg.cpp +++ b/core/src/engine/Hamiltonian_Heisenberg.cpp @@ -27,7 +27,7 @@ namespace Engine intfield anisotropy_indices, scalarfield anisotropy_magnitudes, vectorfield anisotropy_normals, pairfield exchange_pairs, scalarfield exchange_magnitudes, pairfield dmi_pairs, scalarfield dmi_magnitudes, vectorfield dmi_normals, - DDI_Method ddi_method, intfield ddi_n_periodic_images, scalar ddi_radius, + DDI_Method ddi_method, intfield ddi_n_periodic_images, bool ddi_pb_zero_padding, scalar ddi_radius, quadrupletfield quadruplets, scalarfield quadruplet_magnitudes, std::shared_ptr geometry, intfield boundary_conditions @@ -39,7 +39,7 @@ namespace Engine exchange_pairs_in(exchange_pairs), exchange_magnitudes_in(exchange_magnitudes), exchange_shell_magnitudes(0), dmi_pairs_in(dmi_pairs), dmi_magnitudes_in(dmi_magnitudes), dmi_normals_in(dmi_normals), dmi_shell_magnitudes(0), dmi_shell_chirality(0), quadruplets(quadruplets), quadruplet_magnitudes(quadruplet_magnitudes), - ddi_method(ddi_method), ddi_n_periodic_images(ddi_n_periodic_images), ddi_cutoff_radius(ddi_radius), + ddi_method(ddi_method), ddi_n_periodic_images(ddi_n_periodic_images), ddi_pb_zero_padding(ddi_pb_zero_padding) ,ddi_cutoff_radius(ddi_radius), fft_plan_reverse(FFT::FFT_Plan()), fft_plan_spins(FFT::FFT_Plan()) { // Generate interaction pairs, constants etc. @@ -52,7 +52,7 @@ namespace Engine intfield anisotropy_indices, scalarfield anisotropy_magnitudes, vectorfield anisotropy_normals, scalarfield exchange_shell_magnitudes, scalarfield dmi_shell_magnitudes, int dm_chirality, - DDI_Method ddi_method, intfield ddi_n_periodic_images, scalar ddi_radius, + DDI_Method ddi_method, intfield ddi_n_periodic_images, bool ddi_pb_zero_padding, scalar ddi_radius, quadrupletfield quadruplets, scalarfield quadruplet_magnitudes, std::shared_ptr geometry, intfield boundary_conditions @@ -64,7 +64,7 @@ namespace Engine exchange_pairs_in(0), exchange_magnitudes_in(0), exchange_shell_magnitudes(exchange_shell_magnitudes), dmi_pairs_in(0), dmi_magnitudes_in(0), dmi_normals_in(0), dmi_shell_magnitudes(dmi_shell_magnitudes), dmi_shell_chirality(dm_chirality), quadruplets(quadruplets), quadruplet_magnitudes(quadruplet_magnitudes), - ddi_method(ddi_method), ddi_n_periodic_images(ddi_n_periodic_images), ddi_cutoff_radius(ddi_radius), + ddi_method(ddi_method), ddi_n_periodic_images(ddi_n_periodic_images), ddi_pb_zero_padding(ddi_pb_zero_padding), ddi_cutoff_radius(ddi_radius), fft_plan_reverse(FFT::FFT_Plan()), fft_plan_spins(FFT::FFT_Plan()) { @@ -1165,15 +1165,23 @@ namespace Engine if( ddi_method != DDI_Method::FFT ) return; + // We perform zero-padding in a lattice direction if the dimension of the system is greater than 1 *and* + // - the boundary conditions are open, or + // - the boundary conditions are periodic and zero-padding is explicitly requested n_cells_padded.resize(3); - n_cells_padded[0] = (geometry->n_cells[0] > 1) ? 2 * geometry->n_cells[0] : 1; - n_cells_padded[1] = (geometry->n_cells[1] > 1) ? 2 * geometry->n_cells[1] : 1; - n_cells_padded[2] = (geometry->n_cells[2] > 1) ? 2 * geometry->n_cells[2] : 1; + for(int i=0; i<3; i++) + { + n_cells_padded[i] = geometry->n_cells[i]; + bool perform_zero_padding = geometry->n_cells[i] > 1 && (boundary_conditions[i] == 0 || ddi_pb_zero_padding); + if(perform_zero_padding) + n_cells_padded[i] *= 2; + } + sublattice_size = n_cells_padded[0] * n_cells_padded[1] * n_cells_padded[2]; FFT::FFT_Init(); - //workaround for bug in kissfft - //kissfft_ndr does not perform one-dimensional ffts properly + // Workaround for bug in kissfft + // kissfft_ndr does not perform one-dimensional FFTs properly #ifndef SPIRIT_USE_FFTW int number_of_one_dims = 0; for( int i=0; i<3; i++ ) @@ -1185,7 +1193,7 @@ namespace Engine inter_sublattice_lookup.resize(geometry->n_cell_atoms * geometry->n_cell_atoms); - //we dont need to transform over length 1 dims + // We dont need to transform over length 1 dims std::vector fft_dims; for( int i = 2; i >= 0; i-- ) //notice that reverse order is important! { @@ -1193,7 +1201,7 @@ namespace Engine fft_dims.push_back(n_cells_padded[i]); } - //Count how many distinct inter-lattice contributions we need to store + // Count how many distinct inter-lattice contributions we need to store n_inter_sublattice = 0; for( int i = 0; i < geometry->n_cell_atoms; i++ ) { @@ -1204,7 +1212,7 @@ namespace Engine } } - //Create fft plans. + // Create FFT plans FFT::FFT_Plan fft_plan_dipole = FFT::FFT_Plan(fft_dims, false, 6 * n_inter_sublattice, sublattice_size); fft_plan_spins = FFT::FFT_Plan(fft_dims, false, 3 * geometry->n_cell_atoms, sublattice_size); fft_plan_reverse = FFT::FFT_Plan(fft_dims, true, 3 * geometry->n_cell_atoms, sublattice_size); @@ -1230,7 +1238,7 @@ namespace Engine (it_bounds_pointwise_mult[fft_dims.size() - 1] /= 2 )++; #endif - //perform FFT of dipole matrices + // Perform FFT of dipole matrices int img_a = boundary_conditions[0] == 0 ? 0 : ddi_n_periodic_images[0]; int img_b = boundary_conditions[1] == 0 ? 0 : ddi_n_periodic_images[1]; int img_c = boundary_conditions[2] == 0 ? 0 : ddi_n_periodic_images[2]; diff --git a/core/src/engine/Hamiltonian_Heisenberg.cu b/core/src/engine/Hamiltonian_Heisenberg.cu index 40cebb3fe..09b28802f 100644 --- a/core/src/engine/Hamiltonian_Heisenberg.cu +++ b/core/src/engine/Hamiltonian_Heisenberg.cu @@ -28,7 +28,7 @@ namespace Engine intfield anisotropy_indices, scalarfield anisotropy_magnitudes, vectorfield anisotropy_normals, pairfield exchange_pairs, scalarfield exchange_magnitudes, pairfield dmi_pairs, scalarfield dmi_magnitudes, vectorfield dmi_normals, - DDI_Method ddi_method, intfield ddi_n_periodic_images, scalar ddi_radius, + DDI_Method ddi_method, intfield ddi_n_periodic_images, bool ddi_pb_zero_padding, scalar ddi_radius, quadrupletfield quadruplets, scalarfield quadruplet_magnitudes, std::shared_ptr geometry, intfield boundary_conditions @@ -40,7 +40,7 @@ namespace Engine exchange_pairs_in(exchange_pairs), exchange_magnitudes_in(exchange_magnitudes), exchange_shell_magnitudes(0), dmi_pairs_in(dmi_pairs), dmi_magnitudes_in(dmi_magnitudes), dmi_normals_in(dmi_normals), dmi_shell_magnitudes(0), dmi_shell_chirality(0), quadruplets(quadruplets), quadruplet_magnitudes(quadruplet_magnitudes), - ddi_method(ddi_method), ddi_n_periodic_images(ddi_n_periodic_images), ddi_cutoff_radius(ddi_radius), + ddi_method(ddi_method), ddi_n_periodic_images(ddi_n_periodic_images), ddi_pb_zero_padding(ddi_pb_zero_padding), ddi_cutoff_radius(ddi_radius), fft_plan_reverse(FFT::FFT_Plan()), fft_plan_spins(FFT::FFT_Plan()) { // Generate interaction pairs, constants etc. @@ -53,7 +53,7 @@ namespace Engine intfield anisotropy_indices, scalarfield anisotropy_magnitudes, vectorfield anisotropy_normals, scalarfield exchange_shell_magnitudes, scalarfield dmi_shell_magnitudes, int dmi_shell_chirality, - DDI_Method ddi_method, intfield ddi_n_periodic_images, scalar ddi_radius, + DDI_Method ddi_method, intfield ddi_n_periodic_images, bool ddi_pb_zero_padding, scalar ddi_radius, quadrupletfield quadruplets, scalarfield quadruplet_magnitudes, std::shared_ptr geometry, intfield boundary_conditions @@ -65,7 +65,7 @@ namespace Engine exchange_pairs_in(0), exchange_magnitudes_in(0), exchange_shell_magnitudes(exchange_shell_magnitudes), dmi_pairs_in(0), dmi_magnitudes_in(0), dmi_normals_in(0), dmi_shell_magnitudes(dmi_shell_magnitudes), dmi_shell_chirality(dmi_shell_chirality), quadruplets(quadruplets), quadruplet_magnitudes(quadruplet_magnitudes), - ddi_method(ddi_method), ddi_n_periodic_images(ddi_n_periodic_images), ddi_cutoff_radius(ddi_radius), + ddi_method(ddi_method), ddi_n_periodic_images(ddi_n_periodic_images), ddi_pb_zero_padding(ddi_pb_zero_padding), ddi_cutoff_radius(ddi_radius), fft_plan_reverse(FFT::FFT_Plan()), fft_plan_spins(FFT::FFT_Plan()) { // Generate interaction pairs, constants etc. @@ -1038,8 +1038,8 @@ namespace Engine { int tupel[3]; int sublattice_size = iteration_bounds[0] * iteration_bounds[1] * iteration_bounds[2]; - //prefactor of ddi interaction - scalar mult = 2.0133545*1e-28 * 0.057883817555 * 0.057883817555 / ( 4*3.141592653589793238462643383279502884197169399375105820974 * 1e-30 ); + // Prefactor of ddi interaction + scalar mult = C::mu_0 * C::mu_B * C::mu_B / ( 4*C::Pi * 1e-30 ); for(int i = blockIdx.x * blockDim.x + threadIdx.x; i < sublattice_size; i += blockDim.x * gridDim.x) { cu_tupel_from_idx(i, tupel, iteration_bounds, 3); // tupel now is {a, b, c} @@ -1063,7 +1063,7 @@ namespace Engine Vector3 diff; - //iterate over periodic images + // Iterate over periodic images for(int a_pb = - img[0]; a_pb <= img[0]; a_pb++) { for(int b_pb = - img[1]; b_pb <= img[1]; b_pb++) @@ -1144,23 +1144,30 @@ namespace Engine if(ddi_method != DDI_Method::FFT) return; + // We perform zero-padding in a lattice direction if the dimension of the system is greater than 1 *and* + // - the boundary conditions are open, or + // - the boundary conditions are periodic and zero-padding is explicitly requested n_cells_padded.resize(3); - n_cells_padded[0] = (geometry->n_cells[0] > 1) ? 2 * geometry->n_cells[0] : 1; - n_cells_padded[1] = (geometry->n_cells[1] > 1) ? 2 * geometry->n_cells[1] : 1; - n_cells_padded[2] = (geometry->n_cells[2] > 1) ? 2 * geometry->n_cells[2] : 1; + for(int i=0; i<3; i++) + { + n_cells_padded[i] = geometry->n_cells[i]; + bool perform_zero_padding = geometry->n_cells[i] > 1 && (boundary_conditions[i] == 0 || ddi_pb_zero_padding); + if(perform_zero_padding) + n_cells_padded[i] *= 2; + } sublattice_size = n_cells_padded[0] * n_cells_padded[1] * n_cells_padded[2]; inter_sublattice_lookup.resize(geometry->n_cell_atoms * geometry->n_cell_atoms); - //we dont need to transform over length 1 dims + // We dont need to transform over length 1 dims std::vector fft_dims; - for(int i = 2; i >= 0; i--) //notice that reverse order is important! + for(int i = 2; i >= 0; i--) // Notice that reverse order is important! { if(n_cells_padded[i] > 1) fft_dims.push_back(n_cells_padded[i]); } - //Count how many distinct inter-lattice contributions we need to store + // Count how many distinct inter-lattice contributions we need to store n_inter_sublattice = 0; for(int i = 0; i < geometry->n_cell_atoms; i++) { @@ -1171,7 +1178,7 @@ namespace Engine } } - //Set the iteration bounds for the nested for loops that are flattened in the kernels + // Set the iteration bounds for the nested for loops that are flattened in the kernels it_bounds_write_spins = { geometry->n_cell_atoms, geometry->n_cells[0], geometry->n_cells[1], @@ -1200,14 +1207,14 @@ namespace Engine FFT::get_strides(temp_s, {3, this->geometry->n_cell_atoms, n_cells_padded[0], n_cells_padded[1], n_cells_padded[2]}); FFT::get_strides(temp_d, {6, n_inter_sublattice, n_cells_padded[0], n_cells_padded[1], n_cells_padded[2]}); - //perform FFT of dipole matrices + // Perform FFT of dipole matrices int img_a = boundary_conditions[0] == 0 ? 0 : ddi_n_periodic_images[0]; int img_b = boundary_conditions[1] == 0 ? 0 : ddi_n_periodic_images[1]; int img_c = boundary_conditions[2] == 0 ? 0 : ddi_n_periodic_images[2]; FFT_Dipole_Matrices(fft_plan_dipole, img_a, img_b, img_c); transformed_dipole_matrices = std::move(fft_plan_dipole.cpx_ptr); - }//end prepare + }// End prepare void Hamiltonian_Heisenberg::Clean_DDI() { diff --git a/core/src/io/Configparser.cpp b/core/src/io/Configparser.cpp index 11407d0a0..d0e2042fa 100644 --- a/core/src/io/Configparser.cpp +++ b/core/src/io/Configparser.cpp @@ -1135,6 +1135,7 @@ namespace IO auto ddi_method = Engine::DDI_Method::None; intfield ddi_n_periodic_images = { 4, 4, 4 }; scalar ddi_radius = 0.0; + bool ddi_pb_zero_padding = true; // ------------ Quadruplet Interactions ------------ int n_quadruplets = 0; @@ -1347,7 +1348,7 @@ namespace IO // Number of periodical images myfile.Read_3Vector(ddi_n_periodic_images, "ddi_n_periodic_images"); - // myfile.Read_Single(ddi_n_periodic_images, "ddi_n_periodic_images"); + myfile.Read_Single(ddi_pb_zero_padding, "ddi_pb_zero_padding"); // Dipole-dipole cutoff radius myfile.Read_Single(ddi_radius, "ddi_radius"); @@ -1405,6 +1406,7 @@ namespace IO Log(Log_Level::Parameter, Log_Sender::IO, fmt::format(" {:<21} = {}", "ddi_method", ddi_method_str)); Log(Log_Level::Parameter, Log_Sender::IO, fmt::format(" {:<21} = ({} {} {})", "ddi_n_periodic_images", ddi_n_periodic_images[0], ddi_n_periodic_images[1], ddi_n_periodic_images[2])); Log(Log_Level::Parameter, Log_Sender::IO, fmt::format(" {:<21} = {}", "ddi_radius", ddi_radius)); + Log(Log_Level::Parameter, Log_Sender::IO, fmt::format(" {:<21} = {}", "ddi_pb_zero_padding", ddi_pb_zero_padding)); std::unique_ptr hamiltonian; @@ -1415,7 +1417,7 @@ namespace IO anisotropy_index, anisotropy_magnitude, anisotropy_normal, exchange_magnitudes, dmi_magnitudes, dm_chirality, - ddi_method, ddi_n_periodic_images, ddi_radius, + ddi_method, ddi_n_periodic_images, ddi_pb_zero_padding, ddi_radius, quadruplets, quadruplet_magnitudes, geometry, boundary_conditions @@ -1428,7 +1430,7 @@ namespace IO anisotropy_index, anisotropy_magnitude, anisotropy_normal, exchange_pairs, exchange_magnitudes, dmi_pairs, dmi_magnitudes, dmi_normals, - ddi_method, ddi_n_periodic_images, ddi_radius, + ddi_method, ddi_n_periodic_images, ddi_pb_zero_padding, ddi_radius, quadruplets, quadruplet_magnitudes, geometry, boundary_conditions diff --git a/core/test/input/physics_ddi.cfg b/core/test/input/physics_ddi.cfg index e128b57bc..dee8a239d 100644 --- a/core/test/input/physics_ddi.cfg +++ b/core/test/input/physics_ddi.cfg @@ -52,6 +52,7 @@ ddi_n_periodic_images 4 4 4 ### DDI cutoff radius (if cutoff is used) ddi_radius 10.0 +ddi_pb_zero_padding 1 ### Pairs n_interaction_pairs 0 diff --git a/input/input.cfg b/input/input.cfg index 3ca3727e7..2d3c5546c 100644 --- a/input/input.cfg +++ b/input/input.cfg @@ -54,6 +54,9 @@ ddi_n_periodic_images 4 4 4 ### DDI cutoff radius (if cutoff is used) ddi_radius 0.0 +# If set to 1 zero padding is performed even for periodic directions, if set to 0 it is skipped for peridodic directions (improves speed and memory footprint) +ddi_pb_zero_padding 1 + ### Pairs n_interaction_pairs 3 i j da db dc Jij Dij Dijx Dijy Dijz diff --git a/ui-cpp/src/HamiltonianHeisenbergWidget.cpp b/ui-cpp/src/HamiltonianHeisenbergWidget.cpp index 4d03350c4..e4accdcd2 100644 --- a/ui-cpp/src/HamiltonianHeisenbergWidget.cpp +++ b/ui-cpp/src/HamiltonianHeisenbergWidget.cpp @@ -57,6 +57,7 @@ void HamiltonianHeisenbergWidget::Load_Contents() { float d, vd[3], jij[100], dij[100]; int ddi_method, ddi_n_periodic_images[3]; + bool pb_zero_padding; int n_neigh_shells_exchange, n_neigh_shells_dmi, dm_chirality; int n_basis_atoms = Geometry_Get_N_Cell_Atoms(state.get()); std::vector mu_s(n_basis_atoms); @@ -109,7 +110,7 @@ void HamiltonianHeisenbergWidget::Load_Contents() for (int i = 0; i < n_neigh_shells_dmi; ++i) this->dmi_shells[i]->setValue(dij[i]); // DDI - Hamiltonian_Get_DDI(state.get(), &ddi_method, ddi_n_periodic_images, &d); + Hamiltonian_Get_DDI(state.get(), &ddi_method, ddi_n_periodic_images, &d, &pb_zero_padding); this->checkBox_ddi->setChecked( ddi_method != SPIRIT_DDI_METHOD_NONE ); if( ddi_method == SPIRIT_DDI_METHOD_NONE ) this->comboBox_ddi_method->setCurrentIndex(0); @@ -123,6 +124,7 @@ void HamiltonianHeisenbergWidget::Load_Contents() this->spinBox_ddi_n_periodic_b->setValue(ddi_n_periodic_images[1]); this->spinBox_ddi_n_periodic_c->setValue(ddi_n_periodic_images[2]); this->doubleSpinBox_ddi_radius->setValue(d); + this->checkBox_ddi_pb_zero_padding->setChecked(pb_zero_padding); } @@ -464,8 +466,8 @@ void HamiltonianHeisenbergWidget::set_ddi() n_periodic_images[2] = this->spinBox_ddi_n_periodic_c->value(); float radius = this->doubleSpinBox_ddi_radius->value(); - - Hamiltonian_Set_DDI(state.get(), method, n_periodic_images, radius, idx_image); + bool pb_zero_padding = this->checkBox_ddi_pb_zero_padding->isChecked(); + Hamiltonian_Set_DDI(state.get(), method, n_periodic_images, radius, pb_zero_padding, idx_image); }; if (this->comboBox_Hamiltonian_Ani_ApplyTo->currentText() == "Current Image") @@ -556,6 +558,7 @@ void HamiltonianHeisenbergWidget::Setup_Slots() connect(this->spinBox_ddi_n_periodic_b, SIGNAL(editingFinished()), this, SLOT(set_ddi())); connect(this->spinBox_ddi_n_periodic_c, SIGNAL(editingFinished()), this, SLOT(set_ddi())); connect(this->doubleSpinBox_ddi_radius, SIGNAL(editingFinished()), this, SLOT(set_ddi())); + connect(this->checkBox_ddi_pb_zero_padding, SIGNAL(stateChanged(int)), this, SLOT(set_ddi())); // Pairs connect(this->pushButton_pairs_apply, SIGNAL(clicked()), this, SLOT(set_pairs_from_text())); connect(this->pushButton_pairs_fromfile, SIGNAL(clicked()), this, SLOT(set_pairs_from_file())); diff --git a/ui-cpp/ui/HamiltonianHeisenbergWidget.ui b/ui-cpp/ui/HamiltonianHeisenbergWidget.ui index 45b74d60d..c030a9342 100644 --- a/ui-cpp/ui/HamiltonianHeisenbergWidget.ui +++ b/ui-cpp/ui/HamiltonianHeisenbergWidget.ui @@ -537,13 +537,6 @@ - - - - 1000000.000000000000000 - - - @@ -576,6 +569,17 @@ + + + + 1000000.000000000000000 + + + + + + + @@ -624,6 +628,30 @@ + + + + + + Qt::Horizontal + + + + 40 + 20 + + + + + + + + PB zero padding + + + + +