1 // SPDX-FileCopyrightText: Copyright (c) 2017-2024, HONEE contributors. 2 // SPDX-License-Identifier: Apache-2.0 OR BSD-2-Clause 3 4 /// @file 5 /// Time-stepping functions for Navier-Stokes example using PETSc 6 7 #include <ceed.h> 8 #include <petscdmplex.h> 9 #include <petscts.h> 10 11 #include <navierstokes.h> 12 #include "../qfunctions/newtonian_state.h" 13 14 // @brief Insert Boundary values if it's a new time 15 PetscErrorCode UpdateBoundaryValues(User user, Vec Q_loc, PetscReal t) { 16 PetscFunctionBeginUser; 17 if (user->time_bc_set != t) { 18 PetscCall(DMPlexInsertBoundaryValues(user->dm, PETSC_TRUE, Q_loc, t, NULL, NULL, NULL)); 19 user->time_bc_set = t; 20 } 21 PetscFunctionReturn(PETSC_SUCCESS); 22 } 23 24 // RHS (Explicit time-stepper) function setup 25 // This is the RHS of the ODE, given as u_t = G(t,u) 26 // This function takes in a state vector Q and writes into G 27 PetscErrorCode RHS_NS(TS ts, PetscReal t, Vec Q, Vec G, void *user_data) { 28 User user = *(User *)user_data; 29 Ceed ceed = user->ceed; 30 PetscScalar dt; 31 Vec Q_loc = user->Q_loc; 32 PetscMemType q_mem_type; 33 34 PetscFunctionBeginUser; 35 // Update time dependent data 36 PetscCall(UpdateBoundaryValues(user, Q_loc, t)); 37 if (user->phys->solution_time_label) PetscCallCeed(ceed, CeedOperatorSetContextDouble(user->op_rhs_ctx->op, user->phys->solution_time_label, &t)); 38 PetscCall(TSGetTimeStep(ts, &dt)); 39 if (user->phys->timestep_size_label) PetscCallCeed(ceed, CeedOperatorSetContextDouble(user->op_rhs_ctx->op, user->phys->timestep_size_label, &dt)); 40 41 PetscCall(ApplyCeedOperatorGlobalToGlobal(Q, G, user->op_rhs_ctx)); 42 43 PetscCall(VecReadPetscToCeed(Q_loc, &q_mem_type, user->q_ceed)); 44 45 // Inverse of the mass matrix 46 PetscCall(KSPSolve(user->mass_ksp, G, G)); 47 48 PetscCall(VecReadCeedToPetsc(user->q_ceed, q_mem_type, Q_loc)); 49 PetscFunctionReturn(PETSC_SUCCESS); 50 } 51 52 // Surface forces function setup 53 static PetscErrorCode Surface_Forces_NS(DM dm, Vec G_loc, PetscInt num_walls, const PetscInt walls[], PetscScalar *reaction_force) { 54 DMLabel face_label; 55 const PetscScalar *g_array; 56 PetscInt dim = 3; 57 MPI_Comm comm = PetscObjectComm((PetscObject)dm); 58 PetscSection section; 59 60 PetscFunctionBeginUser; 61 PetscCall(DMGetLabel(dm, "Face Sets", &face_label)); 62 PetscCall(VecGetArrayRead(G_loc, &g_array)); 63 for (PetscInt w = 0; w < num_walls; w++) { 64 const PetscInt wall = walls[w], *points; 65 IS wall_is; 66 PetscInt num_points, num_comp = 0; 67 68 PetscCall(DMLabelGetStratumIS(face_label, wall, &wall_is)); 69 if (!wall_is) continue; // No wall points on this process, skip 70 71 PetscCall(DMGetLocalSection(dm, §ion)); 72 PetscCall(PetscSectionGetFieldComponents(section, 0, &num_comp)); 73 PetscCall(ISGetSize(wall_is, &num_points)); 74 PetscCall(ISGetIndices(wall_is, &points)); 75 for (PetscInt i = 0; i < num_points; i++) { 76 const PetscInt p = points[i]; 77 const StateConservative *r; 78 PetscInt dof; 79 80 PetscCall(DMPlexPointLocalRead(dm, p, g_array, &r)); 81 PetscCall(PetscSectionGetDof(section, p, &dof)); 82 for (PetscInt node = 0; node < dof / num_comp; node++) { 83 for (PetscInt j = 0; j < dim; j++) { 84 reaction_force[w * dim + j] -= r[node].momentum[j]; 85 } 86 } 87 } 88 PetscCall(ISRestoreIndices(wall_is, &points)); 89 PetscCall(ISDestroy(&wall_is)); 90 } 91 PetscCallMPI(MPI_Allreduce(MPI_IN_PLACE, reaction_force, dim * num_walls, MPIU_SCALAR, MPI_SUM, comm)); 92 PetscCall(VecRestoreArrayRead(G_loc, &g_array)); 93 PetscFunctionReturn(PETSC_SUCCESS); 94 } 95 96 // Implicit time-stepper function setup 97 PetscErrorCode IFunction_NS(TS ts, PetscReal t, Vec Q, Vec Q_dot, Vec G, void *user_data) { 98 User user = *(User *)user_data; 99 Ceed ceed = user->ceed; 100 PetscScalar dt; 101 Vec Q_loc = user->Q_loc, Q_dot_loc = user->Q_dot_loc, G_loc; 102 PetscMemType q_mem_type, q_dot_mem_type, g_mem_type; 103 104 PetscFunctionBeginUser; 105 // Get local vectors 106 PetscCall(DMGetNamedLocalVector(user->dm, "ResidualLocal", &G_loc)); 107 108 // Update time dependent data 109 PetscCall(UpdateBoundaryValues(user, Q_loc, t)); 110 if (user->phys->solution_time_label) PetscCallCeed(ceed, CeedOperatorSetContextDouble(user->op_ifunction, user->phys->solution_time_label, &t)); 111 PetscCall(TSGetTimeStep(ts, &dt)); 112 if (user->phys->timestep_size_label) PetscCallCeed(ceed, CeedOperatorSetContextDouble(user->op_ifunction, user->phys->timestep_size_label, &dt)); 113 114 // Global-to-local 115 PetscCall(DMGlobalToLocalBegin(user->dm, Q, INSERT_VALUES, Q_loc)); 116 PetscCall(DMGlobalToLocalBegin(user->dm, Q_dot, INSERT_VALUES, Q_dot_loc)); 117 PetscCall(DMGlobalToLocalEnd(user->dm, Q, INSERT_VALUES, Q_loc)); 118 PetscCall(DMGlobalToLocalEnd(user->dm, Q_dot, INSERT_VALUES, Q_dot_loc)); 119 120 // Place PETSc vectors in CEED vectors 121 PetscCall(VecReadPetscToCeed(Q_loc, &q_mem_type, user->q_ceed)); 122 PetscCall(VecReadPetscToCeed(Q_dot_loc, &q_dot_mem_type, user->q_dot_ceed)); 123 PetscCall(VecPetscToCeed(G_loc, &g_mem_type, user->g_ceed)); 124 125 // Apply CEED operator 126 PetscCall(PetscLogEventBegin(FLUIDS_CeedOperatorApply, Q, G, 0, 0)); 127 PetscCall(PetscLogGpuTimeBegin()); 128 PetscCallCeed(user->ceed, CeedOperatorApply(user->op_ifunction, user->q_ceed, user->g_ceed, CEED_REQUEST_IMMEDIATE)); 129 PetscCall(PetscLogGpuTimeEnd()); 130 PetscCall(PetscLogEventEnd(FLUIDS_CeedOperatorApply, Q, G, 0, 0)); 131 132 // Restore vectors 133 PetscCall(VecReadCeedToPetsc(user->q_ceed, q_mem_type, Q_loc)); 134 PetscCall(VecReadCeedToPetsc(user->q_dot_ceed, q_dot_mem_type, Q_dot_loc)); 135 PetscCall(VecCeedToPetsc(user->g_ceed, g_mem_type, G_loc)); 136 137 if (user->app_ctx->sgs_model_type == SGS_MODEL_DATA_DRIVEN) { 138 PetscCall(SgsDDApplyIFunction(user, Q_loc, G_loc)); 139 } 140 141 // Local-to-Global 142 PetscCall(VecZeroEntries(G)); 143 PetscCall(DMLocalToGlobal(user->dm, G_loc, ADD_VALUES, G)); 144 145 // Restore vectors 146 PetscCall(DMRestoreNamedLocalVector(user->dm, "ResidualLocal", &G_loc)); 147 PetscFunctionReturn(PETSC_SUCCESS); 148 } 149 150 PetscErrorCode FormIJacobian_NS(TS ts, PetscReal t, Vec Q, Vec Q_dot, PetscReal shift, Mat J, Mat J_pre, void *user_data) { 151 User user = *(User *)user_data; 152 PetscBool J_is_matceed, J_is_mffd, J_pre_is_matceed, J_pre_is_mffd; 153 154 PetscFunctionBeginUser; 155 PetscCall(PetscObjectTypeCompare((PetscObject)J, MATMFFD, &J_is_mffd)); 156 PetscCall(PetscObjectTypeCompare((PetscObject)J, MATCEED, &J_is_matceed)); 157 PetscCall(PetscObjectTypeCompare((PetscObject)J_pre, MATMFFD, &J_pre_is_mffd)); 158 PetscCall(PetscObjectTypeCompare((PetscObject)J_pre, MATCEED, &J_pre_is_matceed)); 159 160 PetscCall(MatCeedSetContextReal(user->mat_ijacobian, "ijacobian time shift", shift)); 161 162 if (J_is_matceed || J_is_mffd) { 163 PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY)); 164 PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY)); 165 } else PetscCall(MatCeedAssembleCOO(user->mat_ijacobian, J)); 166 167 if (J_pre_is_matceed && J != J_pre) { 168 PetscCall(MatAssemblyBegin(J_pre, MAT_FINAL_ASSEMBLY)); 169 PetscCall(MatAssemblyEnd(J_pre, MAT_FINAL_ASSEMBLY)); 170 } else if (!J_pre_is_matceed && !J_pre_is_mffd && J != J_pre) { 171 PetscCall(MatCeedAssembleCOO(user->mat_ijacobian, J_pre)); 172 } 173 PetscFunctionReturn(PETSC_SUCCESS); 174 } 175 176 PetscErrorCode WriteOutput(User user, Vec Q, PetscInt step_no, PetscScalar time) { 177 Vec Q_loc; 178 char file_path[PETSC_MAX_PATH_LEN]; 179 PetscViewer viewer; 180 181 PetscFunctionBeginUser; 182 if (user->app_ctx->checkpoint_vtk) { 183 // Set up output 184 PetscCall(DMGetLocalVector(user->dm, &Q_loc)); 185 PetscCall(PetscObjectSetName((PetscObject)Q_loc, "StateVec")); 186 PetscCall(VecZeroEntries(Q_loc)); 187 PetscCall(DMGlobalToLocal(user->dm, Q, INSERT_VALUES, Q_loc)); 188 189 // Output 190 PetscCall(PetscSNPrintf(file_path, sizeof file_path, "%s/ns-%03" PetscInt_FMT ".vtu", user->app_ctx->output_dir, step_no)); 191 192 PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)Q), file_path, FILE_MODE_WRITE, &viewer)); 193 PetscCall(VecView(Q_loc, viewer)); 194 PetscCall(PetscViewerDestroy(&viewer)); 195 if (user->dm_viz) { 196 Vec Q_refined, Q_refined_loc; 197 char file_path_refined[PETSC_MAX_PATH_LEN]; 198 PetscViewer viewer_refined; 199 200 PetscCall(DMGetGlobalVector(user->dm_viz, &Q_refined)); 201 PetscCall(DMGetLocalVector(user->dm_viz, &Q_refined_loc)); 202 PetscCall(PetscObjectSetName((PetscObject)Q_refined_loc, "Refined")); 203 204 PetscCall(MatInterpolate(user->interp_viz, Q, Q_refined)); 205 PetscCall(VecZeroEntries(Q_refined_loc)); 206 PetscCall(DMGlobalToLocal(user->dm_viz, Q_refined, INSERT_VALUES, Q_refined_loc)); 207 208 PetscCall( 209 PetscSNPrintf(file_path_refined, sizeof file_path_refined, "%s/nsrefined-%03" PetscInt_FMT ".vtu", user->app_ctx->output_dir, step_no)); 210 211 PetscCall(PetscViewerVTKOpen(PetscObjectComm((PetscObject)Q_refined), file_path_refined, FILE_MODE_WRITE, &viewer_refined)); 212 PetscCall(VecView(Q_refined_loc, viewer_refined)); 213 PetscCall(DMRestoreLocalVector(user->dm_viz, &Q_refined_loc)); 214 PetscCall(DMRestoreGlobalVector(user->dm_viz, &Q_refined)); 215 PetscCall(PetscViewerDestroy(&viewer_refined)); 216 } 217 PetscCall(DMRestoreLocalVector(user->dm, &Q_loc)); 218 } 219 220 // Save data in a binary file for continuation of simulations 221 if (user->app_ctx->add_stepnum2bin) { 222 PetscCall(PetscSNPrintf(file_path, sizeof file_path, "%s/ns-solution-%" PetscInt_FMT ".bin", user->app_ctx->output_dir, step_no)); 223 } else { 224 PetscCall(PetscSNPrintf(file_path, sizeof file_path, "%s/ns-solution.bin", user->app_ctx->output_dir)); 225 } 226 PetscCall(PetscViewerBinaryOpen(user->comm, file_path, FILE_MODE_WRITE, &viewer)); 227 228 time /= user->units->second; // Dimensionalize time back 229 PetscCall(HoneeWriteBinaryVec(viewer, Q, time, step_no)); 230 PetscCall(PetscViewerDestroy(&viewer)); 231 PetscFunctionReturn(PETSC_SUCCESS); 232 } 233 234 // CSV Monitor 235 PetscErrorCode TSMonitor_WallForce(TS ts, PetscInt step_no, PetscReal time, Vec Q, void *ctx) { 236 User user = ctx; 237 Vec G_loc; 238 PetscInt num_wall = user->app_ctx->wall_forces.num_wall, dim = 3; 239 const PetscInt *walls = user->app_ctx->wall_forces.walls; 240 PetscViewer viewer = user->app_ctx->wall_forces.viewer; 241 PetscViewerFormat format = user->app_ctx->wall_forces.viewer_format; 242 PetscScalar *reaction_force; 243 PetscBool is_ascii; 244 245 PetscFunctionBeginUser; 246 if (!viewer) PetscFunctionReturn(PETSC_SUCCESS); 247 PetscCall(DMGetNamedLocalVector(user->dm, "ResidualLocal", &G_loc)); 248 PetscCall(PetscCalloc1(num_wall * dim, &reaction_force)); 249 PetscCall(Surface_Forces_NS(user->dm, G_loc, num_wall, walls, reaction_force)); 250 PetscCall(DMRestoreNamedLocalVector(user->dm, "ResidualLocal", &G_loc)); 251 252 PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &is_ascii)); 253 254 if (is_ascii) { 255 if (format == PETSC_VIEWER_ASCII_CSV && !user->app_ctx->wall_forces.header_written) { 256 PetscCall(PetscViewerASCIIPrintf(viewer, "Step,Time,Wall,ForceX,ForceY,ForceZ\n")); 257 user->app_ctx->wall_forces.header_written = PETSC_TRUE; 258 } 259 for (PetscInt w = 0; w < num_wall; w++) { 260 PetscInt wall = walls[w]; 261 if (format == PETSC_VIEWER_ASCII_CSV) { 262 PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT ",%g,%" PetscInt_FMT ",%g,%g,%g\n", step_no, time, wall, 263 reaction_force[w * dim + 0], reaction_force[w * dim + 1], reaction_force[w * dim + 2])); 264 265 } else { 266 PetscCall(PetscViewerASCIIPrintf(viewer, "Wall %" PetscInt_FMT " Forces: Force_x = %12g, Force_y = %12g, Force_z = %12g\n", wall, 267 reaction_force[w * dim + 0], reaction_force[w * dim + 1], reaction_force[w * dim + 2])); 268 } 269 } 270 } 271 PetscCall(PetscFree(reaction_force)); 272 PetscFunctionReturn(PETSC_SUCCESS); 273 } 274 275 // User provided TS Monitor 276 PetscErrorCode TSMonitor_NS(TS ts, PetscInt step_no, PetscReal time, Vec Q, void *ctx) { 277 User user = ctx; 278 279 PetscFunctionBeginUser; 280 // Print every 'checkpoint_interval' steps 281 if (user->app_ctx->checkpoint_interval <= 0 || step_no % user->app_ctx->checkpoint_interval != 0 || 282 (user->app_ctx->cont_steps == step_no && step_no != 0)) { 283 PetscFunctionReturn(PETSC_SUCCESS); 284 } 285 286 PetscCall(WriteOutput(user, Q, step_no, time)); 287 PetscFunctionReturn(PETSC_SUCCESS); 288 } 289 290 // TS: Create, setup, and solve 291 PetscErrorCode TSSolve_NS(DM dm, User user, AppCtx app_ctx, Physics phys, ProblemData problem, Vec *Q, PetscScalar *f_time, TS *ts) { 292 MPI_Comm comm = user->comm; 293 TSAdapt adapt; 294 PetscScalar final_time; 295 296 PetscFunctionBeginUser; 297 PetscCall(TSCreate(comm, ts)); 298 PetscCall(TSSetDM(*ts, dm)); 299 PetscCall(TSSetApplicationContext(*ts, user)); 300 if (phys->implicit) { 301 PetscCall(TSSetType(*ts, TSBDF)); 302 if (user->op_ifunction) { 303 PetscCall(TSSetIFunction(*ts, NULL, IFunction_NS, &user)); 304 } else { // Implicit integrators can fall back to using an RHSFunction 305 PetscCall(TSSetRHSFunction(*ts, NULL, RHS_NS, &user)); 306 } 307 if (user->mat_ijacobian) { 308 PetscCall(DMTSSetIJacobian(dm, FormIJacobian_NS, &user)); 309 } 310 } else { 311 PetscCheck(user->op_rhs_ctx, comm, PETSC_ERR_ARG_NULL, "Problem does not provide RHSFunction"); 312 PetscCall(TSSetType(*ts, TSRK)); 313 PetscCall(TSRKSetType(*ts, TSRK5F)); 314 PetscCall(TSSetRHSFunction(*ts, NULL, RHS_NS, &user)); 315 } 316 PetscCall(TSSetMaxTime(*ts, 500. * user->units->second)); 317 PetscCall(TSSetExactFinalTime(*ts, TS_EXACTFINALTIME_STEPOVER)); 318 if (app_ctx->test_type == TESTTYPE_NONE) PetscCall(TSSetErrorIfStepFails(*ts, PETSC_FALSE)); 319 PetscCall(TSSetTimeStep(*ts, 1.e-2 * user->units->second)); 320 PetscCall(TSGetAdapt(*ts, &adapt)); 321 PetscCall(TSAdaptSetStepLimits(adapt, 1.e-12 * user->units->second, 1.e2 * user->units->second)); 322 PetscCall(TSSetFromOptions(*ts)); 323 if (user->mat_ijacobian) { 324 if (app_ctx->amat_type && !strcmp(app_ctx->amat_type, MATSHELL)) { 325 SNES snes; 326 KSP ksp; 327 Mat Pmat, Amat; 328 329 PetscCall(TSGetSNES(*ts, &snes)); 330 PetscCall(SNESGetKSP(snes, &ksp)); 331 PetscCall(CreateSolveOperatorsFromMatCeed(ksp, user->mat_ijacobian, PETSC_FALSE, &Amat, &Pmat)); 332 PetscCall(TSSetIJacobian(*ts, user->mat_ijacobian, Pmat, NULL, NULL)); 333 PetscCall(MatDestroy(&Amat)); 334 PetscCall(MatDestroy(&Pmat)); 335 } 336 } 337 user->time_bc_set = -1.0; // require all BCs be updated 338 if (app_ctx->cont_steps) { // continue from previous timestep data 339 PetscInt count; 340 PetscViewer viewer; 341 342 if (app_ctx->cont_time <= 0) { // Legacy files did not include step number and time 343 PetscCall(PetscViewerBinaryOpen(comm, app_ctx->cont_time_file, FILE_MODE_READ, &viewer)); 344 PetscCall(PetscViewerBinaryRead(viewer, &app_ctx->cont_time, 1, &count, PETSC_REAL)); 345 PetscCall(PetscViewerDestroy(&viewer)); 346 PetscCheck(app_ctx->cont_steps != -1, comm, PETSC_ERR_ARG_INCOMP, 347 "-continue step number not specified, but checkpoint file does not contain a step number (likely written by older code version)"); 348 } 349 PetscCall(TSSetTime(*ts, app_ctx->cont_time * user->units->second)); 350 PetscCall(TSSetStepNumber(*ts, app_ctx->cont_steps)); 351 } 352 if (app_ctx->test_type == TESTTYPE_NONE) { 353 PetscCall(TSMonitorSet(*ts, TSMonitor_NS, user, NULL)); 354 } 355 if (app_ctx->wall_forces.viewer) { 356 PetscCall(TSMonitorSet(*ts, TSMonitor_WallForce, user, NULL)); 357 } 358 if (app_ctx->turb_spanstats_enable) { 359 PetscCall(TSMonitorSet(*ts, TSMonitor_TurbulenceStatistics, user, NULL)); 360 CeedScalar previous_time = app_ctx->cont_time * user->units->second; 361 PetscCallCeed(user->ceed, 362 CeedOperatorSetContextDouble(user->spanstats.op_stats_collect_ctx->op, user->spanstats.previous_time_label, &previous_time)); 363 } 364 if (app_ctx->diff_filter_monitor) PetscCall(TSMonitorSet(*ts, TSMonitor_DifferentialFilter, user, NULL)); 365 366 if (app_ctx->sgs_train_enable) { 367 PetscCall(TSMonitorSet(*ts, TSMonitor_SGS_DD_Training, user, NULL)); 368 PetscCall(TSSetPostStep(*ts, TSPostStep_SGS_DD_Training)); 369 } 370 371 if (app_ctx->test_type == TESTTYPE_NONE) PetscCall(PrintRunInfo(user, user->phys, problem, *ts)); 372 // Solve 373 PetscReal start_time; 374 PetscInt start_step; 375 PetscCall(TSGetTime(*ts, &start_time)); 376 PetscCall(TSGetStepNumber(*ts, &start_step)); 377 378 PetscCall(PetscLogDefaultBegin()); // So we can use PetscLogStageGetPerfInfo without -log_view 379 PetscPreLoadBegin(PETSC_FALSE, "Fluids Solve"); 380 PetscCall(TSSetTime(*ts, start_time)); 381 PetscCall(TSSetStepNumber(*ts, start_step)); 382 if (PetscPreLoadingOn) { 383 // LCOV_EXCL_START 384 SNES snes; 385 Vec Q_preload; 386 PetscReal rtol; 387 PetscCall(VecDuplicate(*Q, &Q_preload)); 388 PetscCall(VecCopy(*Q, Q_preload)); 389 PetscCall(TSGetSNES(*ts, &snes)); 390 PetscCall(SNESGetTolerances(snes, NULL, &rtol, NULL, NULL, NULL)); 391 PetscCall(SNESSetTolerances(snes, PETSC_DEFAULT, .99, PETSC_DEFAULT, PETSC_DEFAULT, PETSC_DEFAULT)); 392 PetscCall(TSSetSolution(*ts, Q_preload)); 393 PetscCall(TSStep(*ts)); 394 PetscCall(SNESSetTolerances(snes, PETSC_DEFAULT, rtol, PETSC_DEFAULT, PETSC_DEFAULT, PETSC_DEFAULT)); 395 PetscCall(VecDestroy(&Q_preload)); 396 // LCOV_EXCL_STOP 397 } else { 398 PetscCall(PetscBarrier((PetscObject)*ts)); 399 PetscCall(TSSolve(*ts, *Q)); 400 } 401 PetscPreLoadEnd(); 402 403 PetscCall(TSGetSolveTime(*ts, &final_time)); 404 *f_time = final_time; 405 406 if (app_ctx->test_type == TESTTYPE_NONE) { 407 PetscInt step_no; 408 PetscCall(TSGetStepNumber(*ts, &step_no)); 409 if (user->app_ctx->checkpoint_interval > 0 || user->app_ctx->checkpoint_interval == -1) { 410 PetscCall(WriteOutput(user, *Q, step_no, final_time)); 411 } 412 413 PetscLogStage stage_id; 414 PetscEventPerfInfo stage_perf; 415 416 PetscCall(PetscLogStageGetId("Fluids Solve", &stage_id)); 417 PetscCall(PetscLogStageGetPerfInfo(stage_id, &stage_perf)); 418 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Time taken for solution (sec): %g\n", stage_perf.time)); 419 } 420 PetscFunctionReturn(PETSC_SUCCESS); 421 } 422