diff --git a/driver/level3/level3_thread.c b/driver/level3/level3_thread.c
index aeb5e6ed4e..de29247d4e 100644
--- a/driver/level3/level3_thread.c
+++ b/driver/level3/level3_thread.c
@@ -48,6 +48,10 @@
 #define SWITCH_RATIO 2
 #endif
 
+#ifndef GEMM_PREFERED_SIZE
+#define GEMM_PREFERED_SIZE 1
+#endif
+
 //The array of job_t may overflow the stack.
 //Instead, use malloc to alloc job_t.
 #if MAX_CPU_NUMBER > BLAS3_MEM_ALLOC_THRESHOLD
@@ -510,6 +514,16 @@ static int inner_thread(blas_arg_t *args, BLASLONG *range_m, BLASLONG *range_n,
   return 0;
 }
 
+static int round_up(int remainder, int width, int multiple)
+{
+	if (multiple > remainder || width <= multiple)
+		return width;
+	width = (width + multiple - 1) / multiple;
+	width = width * multiple;
+	return width;
+}
+
+
 static int gemm_driver(blas_arg_t *args, BLASLONG *range_m, BLASLONG
 		       *range_n, FLOAT *sa, FLOAT *sb,
                        BLASLONG nthreads_m, BLASLONG nthreads_n) {
@@ -601,9 +615,14 @@ static int gemm_driver(blas_arg_t *args, BLASLONG *range_m, BLASLONG
   num_parts = 0;
   while (m > 0){
     width = blas_quickdivide(m + nthreads_m - num_parts - 1, nthreads_m - num_parts);
+
+    width = round_up(m, width, GEMM_PREFERED_SIZE);
+
     m -= width;
+
     if (m < 0) width = width + m;
     range_M[num_parts + 1] = range_M[num_parts] + width;
+
     num_parts ++;
   }
   for (i = num_parts; i < MAX_CPU_NUMBER; i++) {
@@ -645,9 +664,12 @@ static int gemm_driver(blas_arg_t *args, BLASLONG *range_m, BLASLONG
       if (width < SWITCH_RATIO) {
         width = SWITCH_RATIO;
       }
+      width = round_up(n, width, GEMM_PREFERED_SIZE);
+
       n -= width;
       if (n < 0) width = width + n;
       range_N[num_parts + 1] = range_N[num_parts] + width;
+
       num_parts ++;
     }
     for (j = num_parts; j < MAX_CPU_NUMBER; j++) {
diff --git a/kernel/x86_64/dgemm_beta_skylakex.c b/kernel/x86_64/dgemm_beta_skylakex.c
index 384e9f60ba..6a824c9b50 100644
--- a/kernel/x86_64/dgemm_beta_skylakex.c
+++ b/kernel/x86_64/dgemm_beta_skylakex.c
@@ -55,6 +55,8 @@ int CNAME(BLASLONG m, BLASLONG n, BLASLONG dummy1, FLOAT beta,
 	return 0;
   }
 
+  if (m == 0 || n == 0)
+	return 0;
 
   c_offset = c;
 
@@ -69,7 +71,7 @@ int CNAME(BLASLONG m, BLASLONG n, BLASLONG dummy1, FLOAT beta,
 
       i = m;
 
-      while (i > 32) {
+      while (i >= 32) {
 	  _mm512_storeu_pd(c_offset1, z_zero);
 	  _mm512_storeu_pd(c_offset1 + 8, z_zero);
 	  _mm512_storeu_pd(c_offset1 + 16, z_zero);
@@ -77,7 +79,7 @@ int CNAME(BLASLONG m, BLASLONG n, BLASLONG dummy1, FLOAT beta,
 	  c_offset1 += 32;
 	  i -= 32;
       }
-      while (i > 8) {
+      while (i >= 8) {
 	  _mm512_storeu_pd(c_offset1, z_zero);
 	  c_offset1 += 8;
 	  i -= 8;
diff --git a/kernel/x86_64/sgemm_beta_skylakex.c b/kernel/x86_64/sgemm_beta_skylakex.c
index 54f9664e9d..4e40acadf8 100644
--- a/kernel/x86_64/sgemm_beta_skylakex.c
+++ b/kernel/x86_64/sgemm_beta_skylakex.c
@@ -55,6 +55,8 @@ int CNAME(BLASLONG m, BLASLONG n, BLASLONG dummy1, FLOAT beta,
 	return 0;
   }
 
+  if (n == 0 || m == 0)
+	return;
 
   c_offset = c;
 
@@ -71,13 +73,13 @@ int CNAME(BLASLONG m, BLASLONG n, BLASLONG dummy1, FLOAT beta,
 
       i = m;
 
-      while (i > 32) {
+      while (i >= 32) {
 	  _mm512_storeu_ps(c_offset1, z_zero);
 	  _mm512_storeu_ps(c_offset1 + 16, z_zero);
 	  c_offset1 += 32;
 	  i -= 32;
       }
-      while (i > 8) {
+      while (i >= 8) {
 	  _mm256_storeu_ps(c_offset1, y_zero);
 	  c_offset1 += 8;
 	  i -= 8;
diff --git a/param.h b/param.h
index e4ec1b2b53..d1b2115845 100644
--- a/param.h
+++ b/param.h
@@ -1627,6 +1627,7 @@ USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 #define SYMV_P  8
 
 #define SWITCH_RATIO	32
+#define GEMM_PREFERED_SIZE	32
 
 #ifdef ARCH_X86