1 // SPDX-FileCopyrightText: Copyright (c) 2017-2024, HONEE contributors. 2 // SPDX-License-Identifier: Apache-2.0 OR BSD-2-Clause 3 4 /// @file 5 /// Miscellaneous utility functions 6 7 #include <ceed.h> 8 #include <petscdm.h> 9 #include <petscsf.h> 10 #include <petscts.h> 11 12 #include <navierstokes.h> 13 #include "../qfunctions/mass.h" 14 15 PetscErrorCode ICs_FixMultiplicity(DM dm, CeedData ceed_data, Honee honee, Vec Q_loc, Vec Q, CeedScalar time) { 16 Ceed ceed = honee->ceed; 17 CeedVector mult_vec; 18 PetscMemType m_mem_type; 19 Vec Multiplicity, Multiplicity_loc; 20 21 PetscFunctionBeginUser; 22 if (honee->phys->ics_time_label) PetscCallCeed(ceed, CeedOperatorSetContextDouble(ceed_data->op_ics_ctx->op, honee->phys->ics_time_label, &time)); 23 PetscCall(ApplyCeedOperatorLocalToGlobal(NULL, Q, ceed_data->op_ics_ctx)); 24 25 PetscCallCeed(ceed, CeedElemRestrictionCreateVector(ceed_data->elem_restr_q, &mult_vec, NULL)); 26 27 // -- Get multiplicity 28 PetscCall(DMGetLocalVector(dm, &Multiplicity_loc)); 29 PetscCall(VecPetscToCeed(Multiplicity_loc, &m_mem_type, mult_vec)); 30 PetscCallCeed(ceed, CeedElemRestrictionGetMultiplicity(ceed_data->elem_restr_q, mult_vec)); 31 PetscCall(VecCeedToPetsc(mult_vec, m_mem_type, Multiplicity_loc)); 32 33 PetscCall(DMGetGlobalVector(dm, &Multiplicity)); 34 PetscCall(VecZeroEntries(Multiplicity)); 35 PetscCall(DMLocalToGlobal(dm, Multiplicity_loc, ADD_VALUES, Multiplicity)); 36 37 // -- Fix multiplicity 38 PetscCall(VecPointwiseDivide(Q, Q, Multiplicity)); 39 PetscCall(VecPointwiseDivide(Q_loc, Q_loc, Multiplicity_loc)); 40 41 PetscCall(DMRestoreLocalVector(dm, &Multiplicity_loc)); 42 PetscCall(DMRestoreGlobalVector(dm, &Multiplicity)); 43 PetscCallCeed(ceed, CeedVectorDestroy(&mult_vec)); 44 PetscFunctionReturn(PETSC_SUCCESS); 45 } 46 47 // Record boundary values from initial condition 48 PetscErrorCode SetBCsFromICs(DM dm, Vec Q, Vec Q_loc) { 49 PetscFunctionBeginUser; 50 { // Capture initial condition values in Qbc 51 Vec Qbc; 52 53 PetscCall(DMGetNamedLocalVector(dm, "Qbc", &Qbc)); 54 PetscCall(VecCopy(Q_loc, Qbc)); 55 PetscCall(VecZeroEntries(Q_loc)); 56 PetscCall(DMGlobalToLocal(dm, Q, INSERT_VALUES, Q_loc)); 57 PetscCall(VecAXPY(Qbc, -1., Q_loc)); 58 PetscCall(DMRestoreNamedLocalVector(dm, "Qbc", &Qbc)); 59 } 60 PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMPlexInsertBoundaryValues_C", DMPlexInsertBoundaryValues_FromICs)); 61 62 { // Set boundary mask to zero out essential BCs 63 Vec boundary_mask, ones; 64 65 PetscCall(DMGetNamedLocalVector(dm, "boundary mask", &boundary_mask)); 66 PetscCall(DMGetGlobalVector(dm, &ones)); 67 PetscCall(VecZeroEntries(boundary_mask)); 68 PetscCall(VecSet(ones, 1.0)); 69 PetscCall(DMGlobalToLocal(dm, ones, INSERT_VALUES, boundary_mask)); 70 PetscCall(DMRestoreNamedLocalVector(dm, "boundary mask", &boundary_mask)); 71 PetscCall(DMRestoreGlobalVector(dm, &ones)); 72 } 73 PetscFunctionReturn(PETSC_SUCCESS); 74 } 75 76 PetscErrorCode DMPlexInsertBoundaryValues_FromICs(DM dm, PetscBool insert_essential, Vec Q_loc, PetscReal time, Vec face_geom_FVM, Vec cell_geom_FVM, 77 Vec grad_FVM) { 78 Vec Qbc, boundary_mask; 79 80 PetscFunctionBeginUser; 81 // Mask (zero) Strong BC entries 82 PetscCall(DMGetNamedLocalVector(dm, "boundary mask", &boundary_mask)); 83 PetscCall(VecPointwiseMult(Q_loc, Q_loc, boundary_mask)); 84 PetscCall(DMRestoreNamedLocalVector(dm, "boundary mask", &boundary_mask)); 85 86 PetscCall(DMGetNamedLocalVector(dm, "Qbc", &Qbc)); 87 PetscCall(VecAXPY(Q_loc, 1., Qbc)); 88 PetscCall(DMRestoreNamedLocalVector(dm, "Qbc", &Qbc)); 89 PetscFunctionReturn(PETSC_SUCCESS); 90 } 91 92 // Compare reference solution values with current test run for CI 93 PetscErrorCode RegressionTest(AppCtx app_ctx, Vec Q) { 94 Vec Qref; 95 PetscViewer viewer; 96 PetscReal error, Qrefnorm; 97 MPI_Comm comm = PetscObjectComm((PetscObject)Q); 98 99 PetscFunctionBeginUser; 100 // Read reference file 101 PetscCall(VecDuplicate(Q, &Qref)); 102 PetscCheck(strcmp(app_ctx->test_file_path, "") != 0, comm, PETSC_ERR_FILE_READ, "File for regression test not given"); 103 PetscCall(PetscViewerBinaryOpen(comm, app_ctx->test_file_path, FILE_MODE_READ, &viewer)); 104 PetscCall(HoneeLoadBinaryVec(viewer, Qref, NULL, NULL)); 105 106 // Compute error with respect to reference solution 107 PetscCall(VecAXPY(Q, -1.0, Qref)); 108 PetscCall(VecNorm(Qref, NORM_MAX, &Qrefnorm)); 109 PetscCall(VecScale(Q, 1. / Qrefnorm)); 110 PetscCall(VecNorm(Q, NORM_MAX, &error)); 111 112 // Check error 113 if (error > app_ctx->test_tol) { 114 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Test failed with error norm %g\n", (double)error)); 115 } 116 117 // Cleanup 118 PetscCall(PetscViewerDestroy(&viewer)); 119 PetscCall(VecDestroy(&Qref)); 120 PetscFunctionReturn(PETSC_SUCCESS); 121 } 122 123 // Get error for problems with exact solutions 124 PetscErrorCode PrintError(CeedData ceed_data, DM dm, Honee honee, Vec Q, PetscScalar final_time) { 125 PetscInt loc_nodes; 126 Vec Q_exact, Q_exact_loc; 127 PetscReal rel_error, norm_error, norm_exact; 128 129 PetscFunctionBeginUser; 130 // Get exact solution at final time 131 PetscCall(DMGetGlobalVector(dm, &Q_exact)); 132 PetscCall(DMGetLocalVector(dm, &Q_exact_loc)); 133 PetscCall(VecGetSize(Q_exact_loc, &loc_nodes)); 134 PetscCall(ICs_FixMultiplicity(dm, ceed_data, honee, Q_exact_loc, Q_exact, final_time)); 135 136 // Get |exact solution - obtained solution| 137 PetscCall(VecNorm(Q_exact, NORM_1, &norm_exact)); 138 PetscCall(VecAXPY(Q, -1.0, Q_exact)); 139 PetscCall(VecNorm(Q, NORM_1, &norm_error)); 140 141 rel_error = norm_error / norm_exact; 142 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Relative Error: %g\n", (double)rel_error)); 143 PetscCall(DMRestoreLocalVector(dm, &Q_exact_loc)); 144 PetscCall(DMRestoreGlobalVector(dm, &Q_exact)); 145 PetscFunctionReturn(PETSC_SUCCESS); 146 } 147 148 // Post-processing 149 PetscErrorCode PostProcess(TS ts, CeedData ceed_data, DM dm, ProblemData problem, Honee honee, Vec Q, PetscScalar final_time) { 150 PetscInt steps; 151 TSConvergedReason reason; 152 153 PetscFunctionBeginUser; 154 // Print relative error 155 if (problem->compute_exact_solution_error && honee->app_ctx->test_type == TESTTYPE_NONE) { 156 PetscCall(PrintError(ceed_data, dm, honee, Q, final_time)); 157 } 158 159 // Print final time and number of steps 160 PetscCall(TSGetStepNumber(ts, &steps)); 161 PetscCall(TSGetConvergedReason(ts, &reason)); 162 if (honee->app_ctx->test_type == TESTTYPE_NONE) { 163 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Time integrator %s on time step %" PetscInt_FMT " with final time %g\n", TSConvergedReasons[reason], 164 steps, (double)final_time)); 165 } 166 167 // Output numerical values from command line 168 PetscCall(VecViewFromOptions(Q, NULL, "-vec_view")); 169 170 // Compare reference solution values with current test run for CI 171 if (honee->app_ctx->test_type == TESTTYPE_SOLVER) { 172 PetscCall(RegressionTest(honee->app_ctx, Q)); 173 } 174 PetscFunctionReturn(PETSC_SUCCESS); 175 } 176 177 // Free a plain data context that was allocated using PETSc; returning libCEED error codes 178 int FreeContextPetsc(void *data) { 179 if (PetscFree(data)) return CeedError(NULL, CEED_ERROR_ACCESS, "PetscFree failed"); 180 return CEED_ERROR_SUCCESS; 181 } 182 183 // @brief Return mass qfunction specification for number of components N 184 PetscErrorCode CreateMassQFunction(Ceed ceed, CeedInt N, CeedInt q_data_size, CeedQFunction *qf) { 185 PetscFunctionBeginUser; 186 switch (N) { 187 case 1: 188 PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, Mass_1, Mass_1_loc, qf)); 189 break; 190 case 2: 191 PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 2, Mass_2, Mass_2_loc, qf)); 192 break; 193 case 3: 194 PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 3, Mass_3, Mass_3_loc, qf)); 195 break; 196 case 4: 197 PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 4, Mass_4, Mass_4_loc, qf)); 198 break; 199 case 5: 200 PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, Mass_5, Mass_5_loc, qf)); 201 break; 202 case 7: 203 PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, Mass_7, Mass_7_loc, qf)); 204 break; 205 case 9: 206 PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, Mass_9, Mass_9_loc, qf)); 207 break; 208 case 12: 209 PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, Mass_12, Mass_12_loc, qf)); 210 break; 211 case 22: 212 PetscCallCeed(ceed, CeedQFunctionCreateInterior(ceed, 1, Mass_22, Mass_22_loc, qf)); 213 break; 214 default: 215 SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, 216 "Mass qfunction of size %" CeedInt_FMT " does not exist.\nA new mass qfunction can be easily added; see source code for pattern.", N); 217 } 218 219 PetscCallCeed(ceed, CeedQFunctionAddInput(*qf, "u", N, CEED_EVAL_INTERP)); 220 PetscCallCeed(ceed, CeedQFunctionAddInput(*qf, "qdata", q_data_size, CEED_EVAL_NONE)); 221 PetscCallCeed(ceed, CeedQFunctionAddOutput(*qf, "v", N, CEED_EVAL_INTERP)); 222 PetscCallCeed(ceed, CeedQFunctionSetUserFlopsEstimate(*qf, N)); 223 PetscFunctionReturn(PETSC_SUCCESS); 224 } 225 226 PetscErrorCode NodalProjectionDataDestroy(NodalProjectionData context) { 227 PetscFunctionBeginUser; 228 if (context == NULL) PetscFunctionReturn(PETSC_SUCCESS); 229 230 PetscCall(DMDestroy(&context->dm)); 231 PetscCall(KSPDestroy(&context->ksp)); 232 233 PetscCall(OperatorApplyContextDestroy(context->l2_rhs_ctx)); 234 235 PetscCall(PetscFree(context)); 236 PetscFunctionReturn(PETSC_SUCCESS); 237 } 238 239 /** 240 @brief Sets the value of an option if it is not already set 241 242 @param[in,out] options `PetscOptions` database, or `NULL` for default global database 243 @param[in] name Name of the option (with `-` prepended to it) 244 @param[in] value Value to set the option to 245 **/ 246 PetscErrorCode HoneeOptionsSetValueDefault(PetscOptions options, const char name[], const char value[]) { 247 PetscBool has_option; 248 249 PetscFunctionBeginUser; 250 PetscCall(PetscOptionsHasName(options, NULL, name, &has_option)); 251 if (!has_option) PetscCall(PetscOptionsSetValue(options, name, value)); 252 PetscFunctionReturn(PETSC_SUCCESS); 253 } 254 255 // Print information about the given simulation run 256 PetscErrorCode PrintRunInfo(Honee honee, Physics phys_ctx, ProblemData problem, TS ts) { 257 Ceed ceed = honee->ceed; 258 MPI_Comm comm = PetscObjectComm((PetscObject)ts); 259 260 PetscFunctionBeginUser; 261 // Header and rank 262 char host_name[PETSC_MAX_PATH_LEN]; 263 PetscMPIInt rank, comm_size; 264 PetscCall(PetscGetHostName(host_name, sizeof host_name)); 265 PetscCallMPI(MPI_Comm_rank(comm, &rank)); 266 PetscCallMPI(MPI_Comm_size(comm, &comm_size)); 267 PetscCall(PetscPrintf(comm, 268 "\n-- Navier-Stokes solver - libCEED + PETSc --\n" 269 " MPI:\n" 270 " Host Name : %s\n" 271 " Total ranks : %d\n", 272 host_name, comm_size)); 273 274 // Problem specific info 275 PetscCall(problem->print_info(honee, problem, honee->app_ctx)); 276 277 // libCEED 278 const char *used_resource; 279 CeedMemType mem_type_backend; 280 PetscCallCeed(ceed, CeedGetResource(honee->ceed, &used_resource)); 281 PetscCallCeed(ceed, CeedGetPreferredMemType(honee->ceed, &mem_type_backend)); 282 PetscCall(PetscPrintf(comm, 283 " libCEED:\n" 284 " libCEED Backend : %s\n" 285 " libCEED Backend MemType : %s\n", 286 used_resource, CeedMemTypes[mem_type_backend])); 287 // PETSc 288 { 289 VecType vec_type; 290 char box_faces_str[PETSC_MAX_PATH_LEN] = "3,3,3"; 291 PetscInt dim; 292 293 PetscCall(DMGetDimension(honee->dm, &dim)); 294 if (dim == 2) box_faces_str[3] = '\0'; 295 PetscCall(PetscOptionsGetString(NULL, NULL, "-dm_plex_box_faces", box_faces_str, sizeof(box_faces_str), NULL)); 296 PetscCall(DMGetVecType(honee->dm, &vec_type)); 297 PetscCall(PetscPrintf(comm, 298 " PETSc:\n" 299 " Box Faces : %s\n" 300 " DM VecType : %s\n" 301 " Time Stepping Scheme : %s\n", 302 box_faces_str, vec_type, phys_ctx->implicit ? "implicit" : "explicit")); 303 } 304 { 305 char pmat_type_str[PETSC_MAX_PATH_LEN]; 306 MatType amat_type, pmat_type; 307 Mat Amat, Pmat; 308 TSIJacobianFn *ijacob_function; 309 310 PetscCall(TSGetIJacobian(ts, &Amat, &Pmat, &ijacob_function, NULL)); 311 PetscCall(MatGetType(Amat, &amat_type)); 312 PetscCall(MatGetType(Pmat, &pmat_type)); 313 314 PetscCall(PetscStrncpy(pmat_type_str, pmat_type, sizeof(pmat_type_str))); 315 if (!strcmp(pmat_type, MATCEED)) { 316 MatType pmat_coo_type; 317 char pmat_coo_type_str[PETSC_MAX_PATH_LEN]; 318 319 PetscCall(MatCeedGetCOOMatType(Pmat, &pmat_coo_type)); 320 PetscCall(PetscSNPrintf(pmat_coo_type_str, sizeof(pmat_coo_type_str), " (COO MatType: %s)", pmat_coo_type)); 321 PetscCall(PetscStrlcat(pmat_type_str, pmat_coo_type_str, sizeof(pmat_type_str))); 322 } 323 if (ijacob_function) { 324 PetscCall(PetscPrintf(comm, 325 " IJacobian A MatType : %s\n" 326 " IJacobian P MatType : %s\n", 327 amat_type, pmat_type_str)); 328 } 329 } 330 if (honee->app_ctx->cont_steps) { 331 PetscCall(PetscPrintf(comm, 332 " Continue:\n" 333 " Filename: : %s\n" 334 " Step: : %" PetscInt_FMT "\n" 335 " Time: : %g\n", 336 honee->app_ctx->cont_file, honee->app_ctx->cont_steps, honee->app_ctx->cont_time)); 337 } 338 // Mesh 339 const PetscInt num_comp_q = 5; 340 PetscInt glob_dofs, owned_dofs, local_dofs; 341 const CeedInt num_P = honee->app_ctx->degree + 1, num_Q = num_P + honee->app_ctx->q_extra; 342 PetscCall(DMGetGlobalVectorInfo(honee->dm, &owned_dofs, &glob_dofs, NULL)); 343 PetscCall(DMGetLocalVectorInfo(honee->dm, &local_dofs, NULL, NULL)); 344 PetscCall(PetscPrintf(comm, 345 " Mesh:\n" 346 " Number of 1D Basis Nodes (P) : %" CeedInt_FMT "\n" 347 " Number of 1D Quadrature Points (Q) : %" CeedInt_FMT "\n" 348 " Global DoFs : %" PetscInt_FMT "\n" 349 " DoFs per node : %" PetscInt_FMT "\n" 350 " Global %" PetscInt_FMT "-DoF nodes : %" PetscInt_FMT "\n", 351 num_P, num_Q, glob_dofs, num_comp_q, num_comp_q, glob_dofs / num_comp_q)); 352 // -- Get Partition Statistics 353 PetscCall(PetscPrintf(comm, " Partition: (min,max,median,max/median)\n")); 354 { 355 PetscInt *gather_buffer = NULL; 356 PetscInt part_owned_dofs[3], part_local_dofs[3], part_boundary_dofs[3], part_neighbors[3]; 357 PetscInt median_index = comm_size % 2 ? comm_size / 2 : comm_size / 2 - 1; 358 if (!rank) PetscCall(PetscMalloc1(comm_size, &gather_buffer)); 359 360 PetscCallMPI(MPI_Gather(&owned_dofs, 1, MPIU_INT, gather_buffer, 1, MPIU_INT, 0, comm)); 361 if (!rank) { 362 PetscCall(PetscSortInt(comm_size, gather_buffer)); 363 part_owned_dofs[0] = gather_buffer[0]; // min 364 part_owned_dofs[1] = gather_buffer[comm_size - 1]; // max 365 part_owned_dofs[2] = gather_buffer[median_index]; // median 366 PetscReal part_owned_dof_ratio = (PetscReal)part_owned_dofs[1] / (PetscReal)part_owned_dofs[2]; 367 PetscCall(PetscPrintf( 368 comm, " Global Vector %" PetscInt_FMT "-DoF nodes : %" PetscInt_FMT ", %" PetscInt_FMT ", %" PetscInt_FMT ", %f\n", num_comp_q, 369 part_owned_dofs[0] / num_comp_q, part_owned_dofs[1] / num_comp_q, part_owned_dofs[2] / num_comp_q, part_owned_dof_ratio)); 370 } 371 372 PetscCallMPI(MPI_Gather(&local_dofs, 1, MPIU_INT, gather_buffer, 1, MPIU_INT, 0, comm)); 373 if (!rank) { 374 PetscCall(PetscSortInt(comm_size, gather_buffer)); 375 part_local_dofs[0] = gather_buffer[0]; // min 376 part_local_dofs[1] = gather_buffer[comm_size - 1]; // max 377 part_local_dofs[2] = gather_buffer[median_index]; // median 378 PetscReal part_local_dof_ratio = (PetscReal)part_local_dofs[1] / (PetscReal)part_local_dofs[2]; 379 PetscCall(PetscPrintf( 380 comm, " Local Vector %" PetscInt_FMT "-DoF nodes : %" PetscInt_FMT ", %" PetscInt_FMT ", %" PetscInt_FMT ", %f\n", num_comp_q, 381 part_local_dofs[0] / num_comp_q, part_local_dofs[1] / num_comp_q, part_local_dofs[2] / num_comp_q, part_local_dof_ratio)); 382 } 383 384 if (comm_size != 1) { 385 PetscInt num_remote_roots_total = 0, num_remote_leaves_total = 0, num_ghost_interface_ranks = 0, num_owned_interface_ranks = 0; 386 { 387 PetscSF sf; 388 PetscMPIInt nrranks, niranks; 389 const PetscInt *roffset, *rmine, *rremote, *ioffset, *irootloc; 390 const PetscMPIInt *rranks, *iranks; 391 PetscCall(DMGetSectionSF(honee->dm, &sf)); 392 PetscCall(PetscSFGetRootRanks(sf, &nrranks, &rranks, &roffset, &rmine, &rremote)); 393 PetscCall(PetscSFGetLeafRanks(sf, &niranks, &iranks, &ioffset, &irootloc)); 394 for (PetscInt i = 0; i < nrranks; i++) { 395 if (rranks[i] == rank) continue; // Ignore same-part global->local transfers 396 num_remote_roots_total += roffset[i + 1] - roffset[i]; 397 num_ghost_interface_ranks++; 398 } 399 for (PetscInt i = 0; i < niranks; i++) { 400 if (iranks[i] == rank) continue; 401 num_remote_leaves_total += ioffset[i + 1] - ioffset[i]; 402 num_owned_interface_ranks++; 403 } 404 } 405 PetscCallMPI(MPI_Gather(&num_remote_roots_total, 1, MPIU_INT, gather_buffer, 1, MPIU_INT, 0, comm)); 406 if (!rank) { 407 PetscCall(PetscSortInt(comm_size, gather_buffer)); 408 part_boundary_dofs[0] = gather_buffer[0]; // min 409 part_boundary_dofs[1] = gather_buffer[comm_size - 1]; // max 410 part_boundary_dofs[2] = gather_buffer[median_index]; // median 411 PetscReal part_shared_dof_ratio = (PetscReal)part_boundary_dofs[1] / (PetscReal)part_boundary_dofs[2]; 412 PetscCall(PetscPrintf( 413 comm, " Ghost Interface %" PetscInt_FMT "-DoF nodes : %" PetscInt_FMT ", %" PetscInt_FMT ", %" PetscInt_FMT ", %f\n", 414 num_comp_q, part_boundary_dofs[0] / num_comp_q, part_boundary_dofs[1] / num_comp_q, part_boundary_dofs[2] / num_comp_q, 415 part_shared_dof_ratio)); 416 } 417 418 PetscCallMPI(MPI_Gather(&num_ghost_interface_ranks, 1, MPIU_INT, gather_buffer, 1, MPIU_INT, 0, comm)); 419 if (!rank) { 420 PetscCall(PetscSortInt(comm_size, gather_buffer)); 421 part_neighbors[0] = gather_buffer[0]; // min 422 part_neighbors[1] = gather_buffer[comm_size - 1]; // max 423 part_neighbors[2] = gather_buffer[median_index]; // median 424 PetscReal part_neighbors_ratio = (PetscReal)part_neighbors[1] / (PetscReal)part_neighbors[2]; 425 PetscCall(PetscPrintf(comm, " Ghost Interface Ranks : %" PetscInt_FMT ", %" PetscInt_FMT ", %" PetscInt_FMT ", %f\n", 426 part_neighbors[0], part_neighbors[1], part_neighbors[2], part_neighbors_ratio)); 427 } 428 429 PetscCallMPI(MPI_Gather(&num_remote_leaves_total, 1, MPIU_INT, gather_buffer, 1, MPIU_INT, 0, comm)); 430 if (!rank) { 431 PetscCall(PetscSortInt(comm_size, gather_buffer)); 432 part_boundary_dofs[0] = gather_buffer[0]; // min 433 part_boundary_dofs[1] = gather_buffer[comm_size - 1]; // max 434 part_boundary_dofs[2] = gather_buffer[median_index]; // median 435 PetscReal part_shared_dof_ratio = (PetscReal)part_boundary_dofs[1] / (PetscReal)part_boundary_dofs[2]; 436 PetscCall(PetscPrintf( 437 comm, " Owned Interface %" PetscInt_FMT "-DoF nodes : %" PetscInt_FMT ", %" PetscInt_FMT ", %" PetscInt_FMT ", %f\n", 438 num_comp_q, part_boundary_dofs[0] / num_comp_q, part_boundary_dofs[1] / num_comp_q, part_boundary_dofs[2] / num_comp_q, 439 part_shared_dof_ratio)); 440 } 441 442 PetscCallMPI(MPI_Gather(&num_owned_interface_ranks, 1, MPIU_INT, gather_buffer, 1, MPIU_INT, 0, comm)); 443 if (!rank) { 444 PetscCall(PetscSortInt(comm_size, gather_buffer)); 445 part_neighbors[0] = gather_buffer[0]; // min 446 part_neighbors[1] = gather_buffer[comm_size - 1]; // max 447 part_neighbors[2] = gather_buffer[median_index]; // median 448 PetscReal part_neighbors_ratio = (PetscReal)part_neighbors[1] / (PetscReal)part_neighbors[2]; 449 PetscCall(PetscPrintf(comm, " Owned Interface Ranks : %" PetscInt_FMT ", %" PetscInt_FMT ", %" PetscInt_FMT ", %f\n", 450 part_neighbors[0], part_neighbors[1], part_neighbors[2], part_neighbors_ratio)); 451 } 452 } 453 454 if (!rank) PetscCall(PetscFree(gather_buffer)); 455 } 456 PetscFunctionReturn(PETSC_SUCCESS); 457 } 458