diff --git a/src/stdlib_intrinsics_dot_product.fypp b/src/stdlib_intrinsics_dot_product.fypp
index ce6188c8a..74bb4b4de 100644
--- a/src/stdlib_intrinsics_dot_product.fypp
+++ b/src/stdlib_intrinsics_dot_product.fypp
@@ -34,10 +34,10 @@ pure module function stdlib_dot_product_${s}$(a,b) result(p)
     n = size(a,kind=ilp)
     r = mod(n,chunk)
 
-    abatch(1:r)       = a(1:r)*${cnjg(t,'b(1:r)')}$
+    abatch(1:r)       = ${cnjg(t,'a(1:r)')}$*b(1:r)
     abatch(r+1:chunk) = zero_${s}$
     do i = r+1, n-r, chunk
-        abatch(1:chunk) = abatch(1:chunk) + a(i:i+chunk-1)*${cnjg(t,'b(i:i+chunk-1)')}$
+        abatch(1:chunk) = abatch(1:chunk) + ${cnjg(t,'a(i:i+chunk-1)')}$*b(i:i+chunk-1)
     end do
 
     p = zero_${s}$
@@ -60,11 +60,11 @@ pure module function stdlib_dot_product_kahan_${s}$(a,b) result(p)
     n = size(a,kind=ilp)
     r = mod(n,chunk)
 
-    abatch(1:r)       = a(1:r)*${cnjg(t,'b(1:r)')}$
+    abatch(1:r)       = ${cnjg(t,'a(1:r)')}$*b(1:r)
     abatch(r+1:chunk) = zero_${s}$
     cbatch = zero_${s}$
     do i = r+1, n-r, chunk
-        call kahan_kernel( a(i:i+chunk-1)*${cnjg(t,'b(i:i+chunk-1)')}$ , abatch(1:chunk) , cbatch(1:chunk) )
+        call kahan_kernel( ${cnjg(t,'a(i:i+chunk-1)')}$*b(i:i+chunk-1) , abatch(1:chunk) , cbatch(1:chunk) )
     end do     
 
     p = zero_${s}$
diff --git a/test/intrinsics/test_intrinsics.fypp b/test/intrinsics/test_intrinsics.fypp
index 8aefe09d3..11cf32fdc 100644
--- a/test/intrinsics/test_intrinsics.fypp
+++ b/test/intrinsics/test_intrinsics.fypp
@@ -246,6 +246,27 @@ subroutine test_dot_product(error)
         call check(error, all(err(:)<tolerance) , "complex dot_product is not accurate" )
         if (allocated(error)) return
     end block
+
+    block ! test for https://github.com/fortran-lang/stdlib/issues/1016
+        ${t}$ :: x(128), y(128)
+        real(${k}$) :: z(128,2)
+        real(${k}$), parameter :: tolerance = epsilon(1._${k}$)*100000
+        real(${k}$) :: err(2)
+        ${t}$ :: p(3)
+
+        call random_number(z)
+        x%re = z(:, 1); x%im = z(:, 2)
+        call random_number(z)
+        y%re = z(:, 1); y%im = z(:, 2)
+        
+        p(1) = dot_product(x,y) ! compiler intrinsic
+        p(2) = stdlib_dot_product_kahan(x,y) ! chunked Kahan dot_product
+        p(3) = stdlib_dot_product(x,y)       ! chunked dot_product
+        err(1:2) = sqrt((p(2:3)%re - p(1)%re)**2 + (p(2:3)%im - p(1)%im)**2)
+        
+        call check(error, all(err(:)<tolerance) , "complex dot_product does not conform to the standard" )
+        if (allocated(error)) return
+    end block
     #:endfor
 
 end subroutine