示例
C Interface:
#include "kml_solver.h" int ierr; int n = 8; int nrhs = 1; // Create matrix A int ia[9] = {0, 2, 4, 6, 7, 8, 10, 12, 14}; int ja[14] = {0, 7, 1, 6, 2, 5, 3, 4, 2, 5, 1, 6, 0, 7}; double a[14] = {1.0, 2.0, -2.0, 3.0, 3.0, 4.0, -4.0, 5.0, 4.0, -6.0, 3.0, 7.0, 2.0, 8.0}; KmlSolverMatrixStore storeA; storeA.indexType = KMLSS_INDEX_INT32; storeA.valueType = KMLSS_VALUE_FP64; storeA.nRow = n; storeA.nCol = n; storeA.format = KMLSS_MATRIX_STORE_CSR; storeA.csr.rowOffset = ia; storeA.csr.colIndex = ja; storeA.csr.value = a; KmlSolverMatrixOption optA; optA.fieldMask = KMLSS_MATRIX_OPTION_TYPE; optA.type = KMLSS_MATRIX_GEN; KmlSolverMatrix *A; ierr = KmlSolverMatrixCreate(&A, &storeA, &optA); if (ierr != KMLSS_NO_ERROR) { printf("ERROR when create A: %d\n", ierr); return 1; } // Create vector b double b[8] = {3.0, 1.0, 7.0, -4.0, 5.0, -2.0, 10.0, 10.0}; KmlSolverMatrixStore storeB; storeB.indexType = KMLSS_INDEX_INT32; storeB.valueType = KMLSS_VALUE_FP64; storeB.nRow = n; storeB.nCol = nrhs; storeB.format = KMLSS_MATRIX_STORE_DENSE_COL_MAJOR; storeB.dense.value = b; storeB.dense.ld = n; KmlSolverMatrixOption optB; optB.fieldMask = KMLSS_MATRIX_OPTION_TYPE; optB.type = KMLSS_MATRIX_GEN; KmlSolverMatrix *B; ierr = KmlSolverMatrixCreate(&B, &storeB, &optB); if (ierr != KMLSS_NO_ERROR) { printf("ERROR when create b: %d\n", ierr); return 1; } // Create vector x double x[8] = {0}; KmlSolverMatrixStore storeX; storeX.indexType = KMLSS_INDEX_INT32; storeX.valueType = KMLSS_VALUE_FP64; storeX.nRow = n; storeX.nCol = nrhs; storeX.format = KMLSS_MATRIX_STORE_DENSE_COL_MAJOR; storeX.dense.value = x; storeX.dense.ld = n; KmlSolverMatrixOption optX; optX.fieldMask = KMLSS_MATRIX_OPTION_TYPE; optX.type = KMLSS_MATRIX_GEN; KmlSolverMatrix *X; ierr = KmlSolverMatrixCreate(&X, &storeX, &optX); if (ierr != KMLSS_NO_ERROR) { printf("ERROR when create x: %d\n", ierr); return 1; } // Init solver KmlDssInitOption opt; opt.fieldMask = KMLDSS_INIT_OPTION_BWR_MODE | KMLDSS_INIT_OPTION_NTHREADS; opt.bwrMode = KMLDSS_BWR_OFF; opt.nThreads = 32; KmlDssSolver *solver; ierr = KmlDssInit(&solver, &opt); if (ierr != KMLSS_NO_ERROR) { printf("ERROR in KmlDssInit: %d\n", ierr); return ierr; } // Analyze KmlDssAnalyzeOption optAnalyze; optAnalyze.fieldMask = KMLDSS_ANALYZE_OPTION_MATCHING_TYPE | KMLDSS_ANALYZE_OPTION_RDR_TYPE | KMLDSS_ANALYZE_OPTION_NTHREADS_RDR; optAnalyze.matchingType = KMLDSS_MATCHING_OFF; optAnalyze.rdrType = KMLDSS_RDR_KRDR; optAnalyze.nThreadsRdr = 8; ierr = KmlDssAnalyze(solver, A, &optAnalyze); if (ierr != KMLSS_NO_ERROR) { printf("ERROR in KmlDssAnalyze: %d\n", ierr); return ierr; } // Factorize KmlDssFactorizeOption optFact; optFact.fieldMask = KMLDSS_FACTORIZE_OPTION_PERTURBATION_THRESHOLD; optFact.perturbationThreshold = 1e-8; ierr = KmlDssFactorize(solver, A, &optFact); if (ierr != KMLSS_NO_ERROR) { printf("ERROR in KmlDssFactorize: %d\n", ierr); return ierr; } // Solve KmlDssSolveOption optSolve; optSolve.fieldMask = KMLDSS_SOLVE_OPTION_SOLVE_STAGE | KMLDSS_SOLVE_OPTION_REFINE_METHOD; optSolve.stage = KMLDSS_SOLVE_ALL; optSolve.refineMethod = KMLDSS_REFINE_OFF; ierr = KmlDssSolve(solver, B, X, &optSolve); if (ierr != KMLSS_NO_ERROR) { printf("ERROR in KmlDssSolve: %d\n", ierr); return ierr; } // Output result x printf("Result of first factorize and solve:\n"); for (int i = 0; i < n; i++) { printf("%lf ", x[i]); } printf("\n"); // Set new values of A double a1[14] = {2.0, 3.0, -3.0, 4.0, 4.0, 5.0, -5.0, 6.0, 5.0, -7.0, 4.0, 8.0, 3.0, 9.0}; KmlSolverMatrixSetValue(A, a1); // Factorize with new values ierr = KmlDssFactorize(solver, A, &optFact); if (ierr != KMLSS_NO_ERROR) { printf("ERROR in KmlDssFactorize: %d\n", ierr); return ierr; } // Set new values of B double b1[8] = {1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0}; KmlSolverMatrixSetValue(B, b1); // Solve with new values ierr = KmlDssSolve(solver, B, X, &optSolve); if (ierr != KMLSS_NO_ERROR) { printf("ERROR in KmlDssSolve: %d\n", ierr); return ierr; } // Output new result x printf("Result of second factorize and solve:\n"); for (int i = 0; i < n; i++) { printf("%lf ", x[i]); } printf("\n"); // Query KmlDssInfo info; info.fieldMask = KMLDSS_INFO_PEAK_MEM; ierr = KmlDssQuery(solver, &info); if (ierr != KMLSS_NO_ERROR) { printf("ERROR in KmlDssQuery: %d\n", ierr); return ierr; } printf("Peak memory is %ld Byte\n", info.peakMem); // Destroy KmlDssClean(&solver); KmlSolverMatrixDestroy(&A); KmlSolverMatrixDestroy(&B); KmlSolverMatrixDestroy(&X);
运行结果:
Result of first factorize and solve: 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 1.000000 Result of second factorize and solve: 0.666667 -0.100000 0.226415 -0.200000 0.166667 0.018868 0.175000 -0.111111 Peak memory is 167775740 Byte
父主题: 求解器函数