diff --git a/src/api.f90 b/src/api.f90
index a53f4fc..f6cc22e 100644
--- a/src/api.f90
+++ b/src/api.f90
@@ -203,6 +203,7 @@ subroutine set_mach(M) bind(c, name='set_mach')
 
         if (M /= MINf) then
             MINf = M
+            MINf1 = M
 
             call comset
             if (MINf>0.0 .and. show_output) write (*, 99003) CPStar, QSTar / QINf
@@ -486,14 +487,14 @@ function get_n_cp() bind(c, name='get_n_cp')
         get_n_cp = N
     end function get_n_cp
 
-    subroutine get_cp(x_out, cp_out, n_points) bind(c, name='get_cp')
+    subroutine get_cp(x_out, y_out, cp_out, n_points) bind(c, name='get_cp')
         use s_xfoil, only: comset
         use m_xoper
         use i_xfoil
         implicit none
 
         integer(c_int), intent(in) :: n_points
-        real(c_float), dimension(n_points), intent(inout) :: x_out, cp_out
+        real(c_float), dimension(n_points), intent(inout) :: x_out, y_out, cp_out
 
         real :: beta, bfac, cpcom, cpinc, den
         integer :: i
@@ -507,6 +508,7 @@ subroutine get_cp(x_out, cp_out, n_points) bind(c, name='get_cp')
             den = beta + bfac * cpinc
             cp_out(i) = cpinc / den
             x_out(i) = X(i)
+            y_out(i) = Y(i)
         enddo
     end subroutine get_cp
 
diff --git a/xfoil/xfoil.py b/xfoil/xfoil.py
index 4334ef6..fb87837 100644
--- a/xfoil/xfoil.py
+++ b/xfoil/xfoil.py
@@ -321,12 +321,15 @@ def get_cp_distribution(self):
         -------
         x : np.array
             X-coordinates
+        y : np.array
+            Y-coordinates
         cp : np.ndarray
             Pressure coefficients at the corresponding x-coordinates
         """
         n = self._lib.get_n_cp()
         x = np.zeros(n, dtype=c_float)
+        y = np.zeros(n, dtype=c_float)
         cp = np.zeros(n, dtype=c_float)
 
-        self._lib.get_cp(x.ctypes.data_as(fptr), cp.ctypes.data_as(fptr), byref(c_int(n)))
-        return x.astype(float), cp.astype(float)
+        self._lib.get_cp(x.ctypes.data_as(fptr), y.ctypes.data_as(fptr), cp.ctypes.data_as(fptr), byref(c_int(n)))
+        return x.astype(float), y.astype(float), cp.astype(float)