                 symmetry_count++;
             }
         }
-
+ 
         //Create fft plans.
         fft_plan_d.dims     = fft_dims;
         fft_plan_d.inverse  = false;

        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_Mats(img_a, img_b, img_c);

    }//end prepare

    // Hamiltonian name as string

        // Quadruplets
    }

    __global__ void CU_Write_FFT_Spin_Input(FFT::FFT_real_type* fft_spin_inputs, const Vector3 * spins, int * iteration_bounds, FFT::StrideContainer spin_stride, scalar * mu_s)
    {
        int nos = iteration_bounds[0] * iteration_bounds[1] * iteration_bounds[2] * iteration_bounds[3];
        int tupel[4];
        int idx_pad;
        for(int idx_orig = blockIdx.x * blockDim.x + threadIdx.x; idx_orig < nos; idx_orig += blockDim.x * gridDim.x)
        {
            tupel_from_idx(idx_orig, tupel, iteration_bounds, 4); //tupel now is {ib, a, b, c}
            idx_pad = tupel[0] * spin_stride.basis + tupel[1] * spin_stride.a + tupel[2] * spin_stride.b + tupel[3] * spin_stride.c;
            fft_spin_inputs[idx_pad                        ] = spins[idx_orig][0] * mu_s[idx_orig];
            fft_spin_inputs[idx_pad + 1 * spin_stride.comp ] = spins[idx_orig][1] * mu_s[idx_orig];
            fft_spin_inputs[idx_pad + 2 * spin_stride.comp ] = spins[idx_orig][2] * mu_s[idx_orig];
        }
    }

    void Hamiltonian_Heisenberg::FFT_spins(const vectorfield & spins)
    {
