134b254c5SRichard Tran Mills static char help[] = "Tests basic creation and destruction of PetscRegressor objects.\n\n"; 234b254c5SRichard Tran Mills 334b254c5SRichard Tran Mills /* 434b254c5SRichard Tran Mills Uses PetscRegressor to train a linear model (that is, linear in its coefficients) 5*789736e1SBarry Smith for a quadratic polynomial data-fitting problem. This is example 3.2 in the first (1996) edition of Michael 634b254c5SRichard Tran Mills T. Heath's "Scientific Computing: An Introductory Survey" textbook. 734b254c5SRichard Tran Mills This example and ex1.c are essentially the same, except the input arrays are mean-centered in ex1.c 834b254c5SRichard Tran Mills and are not in ex2.c. (The data in ex2.c correspond to the data as presented in Heath's example.) 934b254c5SRichard Tran Mills */ 1034b254c5SRichard Tran Mills 1134b254c5SRichard Tran Mills #include <petscregressor.h> 1234b254c5SRichard Tran Mills 1334b254c5SRichard Tran Mills int main(int argc, char **args) 1434b254c5SRichard Tran Mills { 1534b254c5SRichard Tran Mills PetscRegressor regressor; 1634b254c5SRichard Tran Mills PetscMPIInt rank; 1734b254c5SRichard Tran Mills Mat X; 1834b254c5SRichard Tran Mills Vec y, y_predicted, coefficients; 1934b254c5SRichard Tran Mills PetscScalar intercept; 2034b254c5SRichard Tran Mills /* y_array[] and X_array[] are NOT mean-centered; in ex1.c they are! */ 2134b254c5SRichard Tran Mills PetscScalar y_array[5] = {1.0, 0.5, 0, 0.5, 2}; 2234b254c5SRichard Tran Mills PetscScalar X_array[10] = {-1.00000, 1.00000, -0.50000, 0.25000, 0.00000, 0.00000, 0.50000, 0.25000, 1.00000, 1.00000}; 2334b254c5SRichard Tran Mills PetscInt rows_ix[5] = {0, 1, 2, 3, 4}; 2434b254c5SRichard Tran Mills PetscInt cols_ix[2] = {0, 1}; 2534b254c5SRichard Tran Mills 2634b254c5SRichard Tran Mills PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 2734b254c5SRichard Tran Mills PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 2834b254c5SRichard Tran Mills 2934b254c5SRichard Tran Mills PetscCall(VecCreate(PETSC_COMM_WORLD, &y)); 3034b254c5SRichard Tran Mills PetscCall(VecSetSizes(y, PETSC_DECIDE, 5)); 3134b254c5SRichard Tran Mills PetscCall(VecSetFromOptions(y)); 3234b254c5SRichard Tran Mills PetscCall(VecDuplicate(y, &y_predicted)); 3334b254c5SRichard Tran Mills PetscCall(MatCreate(PETSC_COMM_WORLD, &X)); 3434b254c5SRichard Tran Mills PetscCall(MatSetSizes(X, PETSC_DECIDE, PETSC_DECIDE, 5, 2)); 3534b254c5SRichard Tran Mills PetscCall(MatSetFromOptions(X)); 3634b254c5SRichard Tran Mills PetscCall(MatSetUp(X)); 3734b254c5SRichard Tran Mills 3834b254c5SRichard Tran Mills if (!rank) { 3934b254c5SRichard Tran Mills PetscCall(VecSetValues(y, 5, rows_ix, y_array, INSERT_VALUES)); 4034b254c5SRichard Tran Mills PetscCall(MatSetValues(X, 5, rows_ix, 2, cols_ix, X_array, ADD_VALUES)); 4134b254c5SRichard Tran Mills } 4234b254c5SRichard Tran Mills PetscCall(VecAssemblyBegin(y)); 4334b254c5SRichard Tran Mills PetscCall(VecAssemblyEnd(y)); 4434b254c5SRichard Tran Mills PetscCall(MatAssemblyBegin(X, MAT_FINAL_ASSEMBLY)); 4534b254c5SRichard Tran Mills PetscCall(MatAssemblyEnd(X, MAT_FINAL_ASSEMBLY)); 4634b254c5SRichard Tran Mills 4734b254c5SRichard Tran Mills PetscCall(PetscRegressorCreate(PETSC_COMM_WORLD, ®ressor)); 4834b254c5SRichard Tran Mills PetscCall(PetscRegressorSetType(regressor, PETSCREGRESSORLINEAR)); 49*789736e1SBarry Smith PetscCall(PetscRegressorSetFromOptions(regressor)); 5034b254c5SRichard Tran Mills PetscCall(PetscRegressorFit(regressor, X, y)); 5134b254c5SRichard Tran Mills PetscCall(PetscRegressorPredict(regressor, X, y_predicted)); 5234b254c5SRichard Tran Mills PetscCall(PetscRegressorLinearGetIntercept(regressor, &intercept)); 5334b254c5SRichard Tran Mills PetscCall(PetscRegressorLinearGetCoefficients(regressor, &coefficients)); 5434b254c5SRichard Tran Mills 55*789736e1SBarry Smith PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Intercept is %lf\n", (double)intercept)); 5634b254c5SRichard Tran Mills PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Coefficients are\n")); 5734b254c5SRichard Tran Mills PetscCall(VecView(coefficients, PETSC_VIEWER_STDOUT_WORLD)); 5834b254c5SRichard Tran Mills PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Predicted values are\n")); 5934b254c5SRichard Tran Mills PetscCall(VecView(y_predicted, PETSC_VIEWER_STDOUT_WORLD)); 6034b254c5SRichard Tran Mills 6134b254c5SRichard Tran Mills PetscCall(MatDestroy(&X)); 6234b254c5SRichard Tran Mills PetscCall(VecDestroy(&y)); 6334b254c5SRichard Tran Mills PetscCall(VecDestroy(&y_predicted)); 6434b254c5SRichard Tran Mills PetscCall(PetscRegressorDestroy(®ressor)); 6534b254c5SRichard Tran Mills 6634b254c5SRichard Tran Mills PetscCall(PetscFinalize()); 6734b254c5SRichard Tran Mills return 0; 6834b254c5SRichard Tran Mills } 69*789736e1SBarry Smith 70*789736e1SBarry Smith /*TEST 71*789736e1SBarry Smith 72*789736e1SBarry Smith build: 73*789736e1SBarry Smith requires: !complex 74*789736e1SBarry Smith 75*789736e1SBarry Smith test: 76*789736e1SBarry Smith suffix: prefix_tao 77*789736e1SBarry Smith args: -regressor_view 78*789736e1SBarry Smith filter: grep -v "tol: " 79*789736e1SBarry Smith 80*789736e1SBarry Smith test: 81*789736e1SBarry Smith suffix: prefix_ksp 82*789736e1SBarry Smith args: -regressor_view -regressor_linear_use_ksp -regressor_linear_ksp_lsqr_monitor 83*789736e1SBarry Smith 84*789736e1SBarry Smith test: 85*789736e1SBarry Smith requires: suitesparse 86*789736e1SBarry Smith suffix: prefix_ksp_qr 87*789736e1SBarry Smith 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 88*789736e1SBarry Smith TODO: Matrix of type composite does not support checking for transpose 89*789736e1SBarry Smith 90*789736e1SBarry Smith TEST*/ 91