Working on everything but lcvx_bad_scale

Author: govindchari

algebra/cuda/cudss_backend.cu

@@ -586,8 +586,8 @@ static void cudss_solve(LinSysData* linsys_data, QOCOWorkspace* work,
 {
   QOCOInt n = linsys_data->K->n;
   (void)iter_ref_iters;  // No iterative refinement for CUDA backend
-  (void)b;  // b is ignored - RHS is already in d_xyzbuff1 via work->rhs->d_data mapping
-  (void)x;  // x is ignored - solution will be in d_xyzbuff2 via work->xyz->d_data mapping
+  // b is copied to d_xyzbuff1 if provided
+  // x will receive the solution from d_xyzbuff2 after cudssExecute
 
   // During solve phase, ALL data is on GPU - NO CPU-GPU transfers
   // work->rhs->d_data points to d_xyzbuff1 (mapped in map_work_buffers_to_device)
@@ -643,6 +643,22 @@ static void cudss_solve(LinSysData* linsys_data, QOCOWorkspace* work,
                          (status == CUDSS_STATUS_INTERNAL_ERROR) ? "INTERNAL_ERROR" : "UNKNOWN";
     fprintf(stderr, "ERROR: cuDSS solve failed with status %d (%s)\n", (int)status, err_str);
   } else {
+    // Solution is now in d_xyzbuff2 (pointed to by d_xyz_matrix)
+    // Copy solution from d_xyz_matrix (d_xyzbuff2) to x
+    if (x) {
+      cudaPointerAttributes attrs;
+      cudaError_t err = cudaPointerGetAttributes(&attrs, x);
+      if (err == cudaSuccess && attrs.type == cudaMemoryTypeDevice) {
+        // x is on device - device-to-device copy
+        CUDA_CHECK(cudaMemcpy(x, linsys_data->d_xyzbuff2, n * sizeof(QOCOFloat),
+                              cudaMemcpyDeviceToDevice));
+      } else {
+        // x is on host - device-to-host copy
+        CUDA_CHECK(cudaMemcpy(x, linsys_data->d_xyzbuff2, n * sizeof(QOCOFloat),
+                              cudaMemcpyDeviceToHost));
+      }
+    }
+    
     // Debug: print solution after solve
     QOCOFloat* sol_host = (QOCOFloat*)malloc(n * sizeof(QOCOFloat));
     if (sol_host) {
@@ -662,10 +678,22 @@ static void cudss_solve(LinSysData* linsys_data, QOCOWorkspace* work,
   }
 
   // Solution is now in d_xyzbuff2 (pointed to by d_xyz_matrix and work->xyz->d_data)
-  // No need to copy - work->xyz->d_data already points to d_xyzbuff2
+  // Solution has been copied to x (if provided) after successful cudssExecute
   #else
   // cuDSS not available - use fallback: copy solution from RHS (will fail convergence but won't crash)
   CUDA_CHECK(cudaMemcpy(linsys_data->d_xyzbuff2, linsys_data->d_xyzbuff1, n * sizeof(QOCOFloat), cudaMemcpyDeviceToDevice));
+  // Copy to x if provided
+  if (x) {
+    cudaPointerAttributes attrs;
+    cudaError_t err = cudaPointerGetAttributes(&attrs, x);
+    if (err == cudaSuccess && attrs.type == cudaMemoryTypeDevice) {
+      CUDA_CHECK(cudaMemcpy(x, linsys_data->d_xyzbuff2, n * sizeof(QOCOFloat),
+                            cudaMemcpyDeviceToDevice));
+    } else {
+      CUDA_CHECK(cudaMemcpy(x, linsys_data->d_xyzbuff2, n * sizeof(QOCOFloat),
+                            cudaMemcpyDeviceToHost));
+    }
+  }
   #endif
 
   // During solve phase, solution stays on device in d_xyzbuff2 (and work->xyz->d_data)

src/kkt.c

@@ -218,28 +218,177 @@ void initialize_ipm(QOCOSolver* solver)
                                rhs, xyz,
                                solver->settings->iter_ref_iters);
 
-  // Copy x part of solution to x.
+  // Copy x part of solution to x (GPU-to-GPU copy).
+  // During solve phase, get_data_vectorf returns device pointers, so copy_arrayf
+  // will detect both xyz and work->x are on device and use a CUDA kernel.
   copy_arrayf(xyz, get_data_vectorf(solver->work->x), solver->work->data->n);
 
-  // Copy y part of solution to y.
+  // Copy y part of solution to y (GPU-to-GPU copy).
   copy_arrayf(&xyz[solver->work->data->n], get_data_vectorf(solver->work->y),
               solver->work->data->p);
 
-  // Copy z part of solution to z.
+  // Copy z part of solution to z (GPU-to-GPU copy).
   copy_arrayf(&xyz[solver->work->data->n + solver->work->data->p],
               get_data_vectorf(solver->work->z), solver->work->data->m);
 
-  // Copy and negate z part of solution to s.
+  // Copy and negate z part of solution to s (GPU-to-GPU copy).
   copy_and_negate_arrayf(
       &xyz[solver->work->data->n + solver->work->data->p],
       get_data_vectorf(solver->work->s), solver->work->data->m);
 
-      // Note: No sync needed here - get_data_vectorf will return device pointer
-      // during solve phase, so vectors are already on device
-
   // Bring s and z to cone C.
   bring2cone(get_data_vectorf(solver->work->s), solver->work->data);
   bring2cone(get_data_vectorf(solver->work->z), solver->work->data);
+
+  // DEBUG: Print xyz and work->x, work->y, work->z, work->s to verify copies
+  fprintf(stderr, "DEBUG initialize_ipm: Starting debug output\n");
+  fflush(stderr);
+  
+  #ifdef QOCO_ALGEBRA_BACKEND_CUDA
+  fprintf(stderr, "DEBUG initialize_ipm: CUDA backend detected\n");
+  fflush(stderr);
+  {
+    QOCOInt n = solver->work->data->n;
+    QOCOInt p = solver->work->data->p;
+    QOCOInt m = solver->work->data->m;
+    QOCOInt total = n + p + m;
+    
+    fprintf(stderr, "DEBUG initialize_ipm: n=%d, p=%d, m=%d, total=%d\n", n, p, m, total);
+    fflush(stderr);
+    
+    // Allocate host buffers
+    QOCOFloat* xyz_host = (QOCOFloat*)malloc(total * sizeof(QOCOFloat));
+    QOCOFloat* x_host = (QOCOFloat*)malloc(n * sizeof(QOCOFloat));
+    QOCOFloat* y_host = (QOCOFloat*)malloc(p * sizeof(QOCOFloat));
+    QOCOFloat* z_host = (QOCOFloat*)malloc(m * sizeof(QOCOFloat));
+    QOCOFloat* s_host = (QOCOFloat*)malloc(m * sizeof(QOCOFloat));
+    
+    fprintf(stderr, "DEBUG initialize_ipm: Allocated buffers: xyz_host=%p, x_host=%p, y_host=%p, z_host=%p, s_host=%p\n",
+            (void*)xyz_host, (void*)x_host, (void*)y_host, (void*)z_host, (void*)s_host);
+    fflush(stderr);
+    
+    if (xyz_host && x_host && y_host && z_host && s_host) {
+      // Copy from device to host
+      extern void* cudaMemcpy(void* dst, const void* src, size_t count, int kind);
+      #define cudaMemcpyDeviceToHost 2
+      
+      fprintf(stderr, "DEBUG initialize_ipm: Copying from device to host...\n");
+      fflush(stderr);
+      
+      cudaMemcpy(xyz_host, xyz, total * sizeof(QOCOFloat), cudaMemcpyDeviceToHost);
+      cudaMemcpy(x_host, get_data_vectorf(solver->work->x), n * sizeof(QOCOFloat), cudaMemcpyDeviceToHost);
+      cudaMemcpy(y_host, get_data_vectorf(solver->work->y), p * sizeof(QOCOFloat), cudaMemcpyDeviceToHost);
+      cudaMemcpy(z_host, get_data_vectorf(solver->work->z), m * sizeof(QOCOFloat), cudaMemcpyDeviceToHost);
+      cudaMemcpy(s_host, get_data_vectorf(solver->work->s), m * sizeof(QOCOFloat), cudaMemcpyDeviceToHost);
+      
+      fprintf(stderr, "DEBUG initialize_ipm: xyz = {");
+      for (QOCOInt i = 0; i < total; ++i) {
+        fprintf(stderr, "%.17g", xyz_host[i]);
+        if (i != total - 1) fprintf(stderr, ", ");
+      }
+      fprintf(stderr, "}\n");
+      fflush(stderr);
+      
+      fprintf(stderr, "DEBUG initialize_ipm: work->x = {");
+      for (QOCOInt i = 0; i < n; ++i) {
+        fprintf(stderr, "%.17g", x_host[i]);
+        if (i != n - 1) fprintf(stderr, ", ");
+      }
+      fprintf(stderr, "}\n");
+      fflush(stderr);
+      
+      fprintf(stderr, "DEBUG initialize_ipm: work->y = {");
+      for (QOCOInt i = 0; i < p; ++i) {
+        fprintf(stderr, "%.17g", y_host[i]);
+        if (i != p - 1) fprintf(stderr, ", ");
+      }
+      fprintf(stderr, "}\n");
+      fflush(stderr);
+      
+      fprintf(stderr, "DEBUG initialize_ipm: work->z = {");
+      for (QOCOInt i = 0; i < m; ++i) {
+        fprintf(stderr, "%.17g", z_host[i]);
+        if (i != m - 1) fprintf(stderr, ", ");
+      }
+      fprintf(stderr, "}\n");
+      fflush(stderr);
+      
+      fprintf(stderr, "DEBUG initialize_ipm: work->s = {");
+      for (QOCOInt i = 0; i < m; ++i) {
+        fprintf(stderr, "%.17g", s_host[i]);
+        if (i != m - 1) fprintf(stderr, ", ");
+      }
+      fprintf(stderr, "}\n");
+      fflush(stderr);
+      
+      free(xyz_host);
+      free(x_host);
+      free(y_host);
+      free(z_host);
+      free(s_host);
+    } else {
+      fprintf(stderr, "DEBUG initialize_ipm: ERROR - malloc failed!\n");
+      fflush(stderr);
+    }
+  }
+  #else
+  fprintf(stderr, "DEBUG initialize_ipm: Builtin backend (not CUDA)\n");
+  fflush(stderr);
+  {
+    QOCOInt n = solver->work->data->n;
+    QOCOInt p = solver->work->data->p;
+    QOCOInt m = solver->work->data->m;
+    QOCOInt total = n + p + m;
+    QOCOFloat* xyz_ptr = xyz;
+    QOCOFloat* x_ptr = get_data_vectorf(solver->work->x);
+    QOCOFloat* y_ptr = get_data_vectorf(solver->work->y);
+    QOCOFloat* z_ptr = get_data_vectorf(solver->work->z);
+    QOCOFloat* s_ptr = get_data_vectorf(solver->work->s);
+    
+    fprintf(stderr, "DEBUG initialize_ipm: xyz = {");
+    for (QOCOInt i = 0; i < total; ++i) {
+      fprintf(stderr, "%.17g", xyz_ptr[i]);
+      if (i != total - 1) fprintf(stderr, ", ");
+    }
+    fprintf(stderr, "}\n");
+    fflush(stderr);
+    
+    fprintf(stderr, "DEBUG initialize_ipm: work->x = {");
+    for (QOCOInt i = 0; i < n; ++i) {
+      fprintf(stderr, "%.17g", x_ptr[i]);
+      if (i != n - 1) fprintf(stderr, ", ");
+    }
+    fprintf(stderr, "}\n");
+    fflush(stderr);
+    
+    fprintf(stderr, "DEBUG initialize_ipm: work->y = {");
+    for (QOCOInt i = 0; i < p; ++i) {
+      fprintf(stderr, "%.17g", y_ptr[i]);
+      if (i != p - 1) fprintf(stderr, ", ");
+    }
+    fprintf(stderr, "}\n");
+    fflush(stderr);
+    
+    fprintf(stderr, "DEBUG initialize_ipm: work->z = {");
+    for (QOCOInt i = 0; i < m; ++i) {
+      fprintf(stderr, "%.17g", z_ptr[i]);
+      if (i != m - 1) fprintf(stderr, ", ");
+    }
+    fprintf(stderr, "}\n");
+    fflush(stderr);
+    
+    fprintf(stderr, "DEBUG initialize_ipm: work->s = {");
+    for (QOCOInt i = 0; i < m; ++i) {
+      fprintf(stderr, "%.17g", s_ptr[i]);
+      if (i != m - 1) fprintf(stderr, ", ");
+    }
+    fprintf(stderr, "}\n");
+    fflush(stderr);
+  }
+  #endif
+  
+  fprintf(stderr, "DEBUG initialize_ipm: Finished debug output\n");
+  fflush(stderr);
 }
 
 void compute_kkt_residual(QOCOProblemData* data, QOCOFloat* x, QOCOFloat* y,

src/qoco_api.c

@@ -197,7 +197,7 @@ void qoco_set_csc(QOCOCscMatrix* A, QOCOInt m, QOCOInt n, QOCOInt Annz,
 
 void set_default_settings(QOCOSettings* settings)
 {
-  settings->max_iters = 1;
+  settings->max_iters = 200;
   settings->bisect_iters = 5;
   settings->ruiz_iters = 0;
   settings->iter_ref_iters = 1;