1 // SPDX-FileCopyrightText: Copyright (c) 2017-2024, HONEE contributors. 2 // SPDX-License-Identifier: Apache-2.0 OR BSD-2-Clause 3 4 /// @file 5 /// Utility functions for setting up Channel flow 6 7 #include "../qfunctions/channel.h" 8 9 #include <ceed.h> 10 #include <petscdm.h> 11 12 #include <navierstokes.h> 13 14 PetscErrorCode NS_CHANNEL(ProblemData problem, DM dm, void *ctx, SimpleBC bc) { 15 Honee honee = *(Honee *)ctx; 16 MPI_Comm comm = honee->comm; 17 Ceed ceed = honee->ceed; 18 ChannelContext channel_ctx; 19 NewtonianIdealGasContext newtonian_ig_ctx; 20 CeedQFunctionContext channel_qfctx; 21 22 PetscFunctionBeginUser; 23 PetscCall(NS_NEWTONIAN_IG(problem, dm, ctx, bc)); 24 PetscCall(PetscCalloc1(1, &channel_ctx)); 25 26 // ------------------------------------------------------ 27 // SET UP Channel 28 // ------------------------------------------------------ 29 PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->ics.qfctx)); 30 problem->ics.qf_func_ptr = ICsChannel; 31 problem->ics.qf_loc = ICsChannel_loc; 32 if (honee->phys->state_var == STATEVAR_CONSERVATIVE) { 33 problem->apply_inflow.qf_func_ptr = Channel_Inflow; 34 problem->apply_inflow.qf_loc = Channel_Inflow_loc; 35 problem->apply_outflow.qf_func_ptr = Channel_Outflow; 36 problem->apply_outflow.qf_loc = Channel_Outflow_loc; 37 } 38 39 // -- Command Line Options 40 CeedScalar umax = 10.; // m/s 41 CeedScalar theta0 = 300.; // K 42 CeedScalar P0 = 1.e5; // Pa 43 PetscReal body_force_scale = 1.; 44 PetscOptionsBegin(comm, NULL, "Options for CHANNEL problem", NULL); 45 PetscCall(PetscOptionsScalar("-umax", "Centerline velocity of the Channel", NULL, umax, &umax, NULL)); 46 PetscCall(PetscOptionsScalar("-theta0", "Wall temperature", NULL, theta0, &theta0, NULL)); 47 PetscCall(PetscOptionsScalar("-P0", "Pressure at outflow", NULL, P0, &P0, NULL)); 48 PetscCall(PetscOptionsReal("-body_force_scale", "Multiplier for body force", NULL, body_force_scale = 1, &body_force_scale, NULL)); 49 PetscOptionsEnd(); 50 51 PetscScalar meter = honee->units->meter; 52 PetscScalar second = honee->units->second; 53 PetscScalar Kelvin = honee->units->Kelvin; 54 PetscScalar Pascal = honee->units->Pascal; 55 56 theta0 *= Kelvin; 57 P0 *= Pascal; 58 umax *= meter / second; 59 60 //-- Setup Problem information 61 CeedScalar H, center; 62 { 63 PetscReal domain_min[3], domain_max[3], domain_size[3]; 64 PetscCall(DMGetBoundingBox(dm, domain_min, domain_max)); 65 for (PetscInt i = 0; i < 3; i++) domain_size[i] = domain_max[i] - domain_min[i]; 66 67 H = 0.5 * domain_size[1] * meter; 68 center = H + domain_min[1] * meter; 69 } 70 71 // Some properties depend on parameters from NewtonianIdealGas 72 PetscCallCeed(ceed, CeedQFunctionContextGetData(problem->apply_vol_rhs.qfctx, CEED_MEM_HOST, &newtonian_ig_ctx)); 73 74 channel_ctx->center = center; 75 channel_ctx->H = H; 76 channel_ctx->theta0 = theta0; 77 channel_ctx->P0 = P0; 78 channel_ctx->umax = umax; 79 channel_ctx->implicit = honee->phys->implicit; 80 channel_ctx->B = body_force_scale * 2 * umax * newtonian_ig_ctx->mu / (H * H); 81 82 { 83 // Calculate Body force 84 CeedScalar cv = newtonian_ig_ctx->cv, cp = newtonian_ig_ctx->cp; 85 CeedScalar Rd = cp - cv; 86 CeedScalar rho = P0 / (Rd * theta0); 87 CeedScalar g[] = {channel_ctx->B / rho, 0., 0.}; 88 PetscCall(PetscArraycpy(newtonian_ig_ctx->g, g, 3)); 89 } 90 channel_ctx->newtonian_ctx = *newtonian_ig_ctx; 91 PetscCallCeed(ceed, CeedQFunctionContextRestoreData(problem->apply_vol_rhs.qfctx, &newtonian_ig_ctx)); 92 93 PetscCallCeed(ceed, CeedQFunctionContextCreate(honee->ceed, &channel_qfctx)); 94 PetscCallCeed(ceed, CeedQFunctionContextSetData(channel_qfctx, CEED_MEM_HOST, CEED_USE_POINTER, sizeof(*channel_ctx), channel_ctx)); 95 PetscCallCeed(ceed, CeedQFunctionContextSetDataDestroy(channel_qfctx, CEED_MEM_HOST, FreeContextPetsc)); 96 97 problem->ics.qfctx = channel_qfctx; 98 PetscCallCeed(ceed, CeedQFunctionContextReferenceCopy(channel_qfctx, &problem->apply_inflow.qfctx)); 99 PetscCallCeed(ceed, CeedQFunctionContextReferenceCopy(channel_qfctx, &problem->apply_outflow.qfctx)); 100 PetscFunctionReturn(PETSC_SUCCESS); 101 } 102