Actual source code: ex2.c

  1: static char help[] = "Tests basic creation and destruction of PetscRegressor objects.\n\n";

  3: /*
  4:     Uses PetscRegressor to train a linear model (that is, linear in its coefficients)
  5:     for a quadratic polynomial data-fitting problem. This is example 3.2 in the first (1996) edition of Michael
  6:     T. Heath's "Scientific Computing: An Introductory Survey" textbook.
  7:     This example and ex1.c are essentially the same, except the input arrays are mean-centered in ex1.c
  8:     and are not in ex2.c. (The data in ex2.c correspond to the data as presented in Heath's example.)
  9: */

 11: #include <petscregressor.h>

 13: int main(int argc, char **args)
 14: {
 15:   PetscRegressor regressor;
 16:   PetscMPIInt    rank;
 17:   Mat            X;
 18:   Vec            y, y_predicted, coefficients;
 19:   PetscScalar    intercept;
 20:   /* y_array[] and X_array[] are NOT mean-centered; in ex1.c they are! */
 21:   PetscScalar y_array[5]  = {1.0, 0.5, 0, 0.5, 2};
 22:   PetscScalar X_array[10] = {-1.00000, 1.00000, -0.50000, 0.25000, 0.00000, 0.00000, 0.50000, 0.25000, 1.00000, 1.00000};
 23:   PetscInt    rows_ix[5]  = {0, 1, 2, 3, 4};
 24:   PetscInt    cols_ix[2]  = {0, 1};

 26:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
 27:   PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));

 29:   PetscCall(VecCreate(PETSC_COMM_WORLD, &y));
 30:   PetscCall(VecSetSizes(y, PETSC_DECIDE, 5));
 31:   PetscCall(VecSetFromOptions(y));
 32:   PetscCall(VecDuplicate(y, &y_predicted));
 33:   PetscCall(MatCreate(PETSC_COMM_WORLD, &X));
 34:   PetscCall(MatSetSizes(X, PETSC_DECIDE, PETSC_DECIDE, 5, 2));
 35:   PetscCall(MatSetFromOptions(X));
 36:   PetscCall(MatSetUp(X));

 38:   if (!rank) {
 39:     PetscCall(VecSetValues(y, 5, rows_ix, y_array, INSERT_VALUES));
 40:     PetscCall(MatSetValues(X, 5, rows_ix, 2, cols_ix, X_array, ADD_VALUES));
 41:   }
 42:   PetscCall(VecAssemblyBegin(y));
 43:   PetscCall(VecAssemblyEnd(y));
 44:   PetscCall(MatAssemblyBegin(X, MAT_FINAL_ASSEMBLY));
 45:   PetscCall(MatAssemblyEnd(X, MAT_FINAL_ASSEMBLY));

 47:   PetscCall(PetscRegressorCreate(PETSC_COMM_WORLD, &regressor));
 48:   PetscCall(PetscRegressorSetType(regressor, PETSCREGRESSORLINEAR));
 49:   PetscCall(PetscRegressorSetFromOptions(regressor));
 50:   PetscCall(PetscRegressorFit(regressor, X, y));
 51:   PetscCall(PetscRegressorPredict(regressor, X, y_predicted));
 52:   PetscCall(PetscRegressorLinearGetIntercept(regressor, &intercept));
 53:   PetscCall(PetscRegressorLinearGetCoefficients(regressor, &coefficients));

 55:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Intercept is %lf\n", (double)intercept));
 56:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Coefficients are\n"));
 57:   PetscCall(VecView(coefficients, PETSC_VIEWER_STDOUT_WORLD));
 58:   PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Predicted values are\n"));
 59:   PetscCall(VecView(y_predicted, PETSC_VIEWER_STDOUT_WORLD));

 61:   PetscCall(MatDestroy(&X));
 62:   PetscCall(VecDestroy(&y));
 63:   PetscCall(VecDestroy(&y_predicted));
 64:   PetscCall(PetscRegressorDestroy(&regressor));

 66:   PetscCall(PetscFinalize());
 67:   return 0;
 68: }

 70: /*TEST

 72:   build:
 73:     requires: !complex !single !__float128 !defined(PETSC_USE_64BIT_INDICES)

 75:   test:
 76:     suffix: prefix_tao
 77:     args: -regressor_view

 79:   test:
 80:     suffix: prefix_ksp
 81:     args: -regressor_view -regressor_linear_use_ksp -regressor_linear_ksp_lsqr_monitor

 83:   test:
 84:     requires: suitesparse
 85:     suffix: prefix_ksp_qr
 86:     args: -regressor_view -regressor_linear_use_ksp -regressor_linear_ksp_lsqr_monitor -regressor_linear_pc_type qr regressor_linear_pc_factor_mat_solver_type spqr
 87:     TODO: Matrix of type composite does not support checking for transpose

 89: TEST*/