Lines Matching +full:- +full:f
1 // Copyright (c) 2017-2026, Lawrence Livermore National Security, LLC and other CEED contributors.
2 // All Rights Reserved. See the top-level LICENSE and NOTICE files for details.
4 // SPDX-License-Identifier: BSD-2-Clause
23 PetscScalar *r, f[4], h[4]; in CompressibleBlasiusResidual() local
24 PetscInt N = blasius->n_cheb; in CompressibleBlasiusResidual()
25 State S_infty = blasius->S_infty; in CompressibleBlasiusResidual()
29 …PetscScalar Ma = Mach(&blasius->newtonian_ctx, S_infty.Y.temperature, U_infty), Pr = Prandtl(&blas… in CompressibleBlasiusResidual()
30 gamma = HeatCapacityRatio(&blasius->newtonian_ctx); in CompressibleBlasiusResidual()
36 // Left boundary conditions f = f' = 0 in CompressibleBlasiusResidual()
37 ChebyshevEval(N, Tf, -1., blasius->eta_max, f); in CompressibleBlasiusResidual()
38 r[0] = f[0]; in CompressibleBlasiusResidual()
39 r[1] = f[1]; in CompressibleBlasiusResidual()
41 // f - right end boundary condition in CompressibleBlasiusResidual()
42 ChebyshevEval(N, Tf, 1., blasius->eta_max, f); in CompressibleBlasiusResidual()
43 r[2] = f[1] - 1.; in CompressibleBlasiusResidual()
45 for (int i = 0; i < N - 3; i++) { in CompressibleBlasiusResidual()
46 ChebyshevEval(N, Tf, blasius->X[i], blasius->eta_max, f); in CompressibleBlasiusResidual()
47 ChebyshevEval(N - 1, Th, blasius->X[i], blasius->eta_max, h); in CompressibleBlasiusResidual()
53 const PetscScalar rho_tilde[2] = {1 / h[0], -h[1] / PetscSqr(h[0])}; in CompressibleBlasiusResidual()
58 r[3 + i] = 2 * (mu_rho_tilde[0] * f[3] + mu_rho_tilde[1] * f[2]) + f[2] * f[0]; in CompressibleBlasiusResidual()
59 …0] * h[2] + mu_rho_tilde[1] * h[1]) + Pr * f[0] * h[1] + Pr * (gamma - 1) * mu_rho_tilde[0] * Pets… in CompressibleBlasiusResidual()
62 // h - left end boundary condition in CompressibleBlasiusResidual()
63 ChebyshevEval(N - 1, Th, -1., blasius->eta_max, h); in CompressibleBlasiusResidual()
64 r[N] = h[0] - blasius->T_wall / S_infty.Y.temperature; in CompressibleBlasiusResidual()
66 // h - right end boundary condition in CompressibleBlasiusResidual()
67 ChebyshevEval(N - 1, Th, 1., blasius->eta_max, h); in CompressibleBlasiusResidual()
68 r[N + 1] = h[0] - 1.; in CompressibleBlasiusResidual()
80 PetscInt N = blasius->n_cheb; in ComputeChebyshevCoefficients()
86 PetscCall(PetscMalloc2(N - 3, &blasius->X, N - 3, &w)); in ComputeChebyshevCoefficients()
87 PetscCall(PetscDTGaussQuadrature(N - 3, -1., 1., blasius->X, w)); in ComputeChebyshevCoefficients()
92 PetscCall(VecSetSizes(sol, PETSC_DECIDE, 2 * N - 1)); in ComputeChebyshevCoefficients()
106 for (int i = 0; i < N; i++) blasius->Tf_cheb[i] = cheb_coefs[i]; in ComputeChebyshevCoefficients()
107 for (int i = 0; i < N - 1; i++) blasius->Th_cheb[i] = cheb_coefs[i + N]; in ComputeChebyshevCoefficients()
110 PetscCall(PetscFree2(blasius->X, w)); in ComputeChebyshevCoefficients()
177 for (PetscInt i = 0; i < 3; i++) domain_size[i] = domain_max[i] - domain_min[i]; in ModifyMesh()
189 PetscCall(PetscOptionsGetIntArray(NULL, NULL, "-dm_plex_box_faces", faces, &nmax, NULL)); in ModifyMesh()
195 PetscReal dy1 = refine_height * (growth - 1) / (pow(growth, N) - 1); in ModifyMesh()
198 PetscReal logdy = (log(domain_max[1]) - log(refine_height)) / (faces[1] - N); in ModifyMesh()
208 …(1 - (coords[i][0] - domain_min[0]) * angle_coeff / domain_max[1]) * dy1 * (pow(growth, coords[i][… in ModifyMesh()
210 PetscInt j = y_box_index - N; in ModifyMesh()
211 …coords[i][1] = (1 - (coords[i][0] - domain_min[0]) * angle_coeff / domain_max[1]) * exp(log(refine… in ModifyMesh()
230 // Determine which y-node we're at in ModifyMesh()
232 …coords[i][1] = (1 - (coords[i][0] - domain_min[0]) * angle_coeff / max_y) * (*node_locs)[y… in ModifyMesh()
243 MPI_Comm comm = user->comm; in NS_BLASIUS()
244 Ceed ceed = user->ceed; in NS_BLASIUS()
254 // ------------------------------------------------------ in NS_BLASIUS()
256 // ------------------------------------------------------ in NS_BLASIUS()
257 problem->ics.qfunction = ICsBlasius; in NS_BLASIUS()
258 problem->ics.qfunction_loc = ICsBlasius_loc; in NS_BLASIUS()
263 CeedScalar delta0 = 4.2e-3; // m in NS_BLASIUS()
267 PetscReal mesh_refine_height = 5.9e-4; // m in NS_BLASIUS()
268 PetscReal mesh_growth = 1.08; // [-] in NS_BLASIUS()
269 PetscInt mesh_Ndelta = 45; // [-] in NS_BLASIUS()
275 …PetscCall(PetscOptionsBool("-weakT", "Change from rho weak to T weak at inflow", NULL, weakT, &wea… in NS_BLASIUS()
276 …PetscCall(PetscOptionsScalar("-velocity_infinity", "Velocity at boundary layer edge", NULL, U_inf,… in NS_BLASIUS()
277 …PetscCall(PetscOptionsScalar("-temperature_infinity", "Temperature at boundary layer edge", NULL, … in NS_BLASIUS()
278 …PetscCall(PetscOptionsHasName(NULL, NULL, "-P0", &P0_set)); // For maintaining behavior of -P0 fl… in NS_BLASIUS()
280 PetscOptionsDeprecated("-P0", "-pressure_infinity", "libCEED 0.12.0", in NS_BLASIUS()
281 …"Use -pressure_infinity to set pressure at boundary layer edge and -idl_pressure to set the IDL re… in NS_BLASIUS()
282 …PetscCall(PetscOptionsScalar("-pressure_infinity", "Pressure at boundary layer edge", NULL, P_inf,… in NS_BLASIUS()
283 …PetscCall(PetscOptionsScalar("-temperature_wall", "Temperature at wall", NULL, T_wall, &T_wall, NU… in NS_BLASIUS()
284 …PetscCall(PetscOptionsScalar("-delta0", "Boundary layer height at inflow", NULL, delta0, &delta0, … in NS_BLASIUS()
285 PetscCall(PetscOptionsInt("-n_chebyshev", "Number of Chebyshev terms", NULL, N, &N, NULL)); in NS_BLASIUS()
286 …PetscCheck(3 <= N && N <= BLASIUS_MAX_N_CHEBYSHEV, comm, PETSC_ERR_ARG_OUTOFRANGE, "-n_chebyshev %… in NS_BLASIUS()
288 if (user->app_ctx->mesh_transform_type == MESH_TRANSFORM_PLATEMESH) { in NS_BLASIUS()
289 …PetscCall(PetscOptionsBoundedInt("-platemesh_Ndelta", "Velocity at boundary layer edge", NULL, mes… in NS_BLASIUS()
290 …PetscCall(PetscOptionsScalar("-platemesh_refine_height", "Height of boundary layer mesh refinement… in NS_BLASIUS()
292 …PetscCall(PetscOptionsScalar("-platemesh_growth", "Geometric growth rate of boundary layer mesh", … in NS_BLASIUS()
293 …PetscCall(PetscOptionsScalar("-platemesh_top_angle", "Geometric top_angle rate of boundary layer m… in NS_BLASIUS()
295 PetscCall(PetscOptionsString("-platemesh_y_node_locs_path", in NS_BLASIUS()
300 …PetscCall(PetscOptionsBool("-stg_use", "Use STG inflow boundary condition", NULL, use_stg, &use_st… in NS_BLASIUS()
303 PetscScalar meter = user->units->meter; in NS_BLASIUS()
304 PetscScalar second = user->units->second; in NS_BLASIUS()
305 PetscScalar Kelvin = user->units->Kelvin; in NS_BLASIUS()
306 PetscScalar Pascal = user->units->Pascal; in NS_BLASIUS()
314 if (user->app_ctx->mesh_transform_type == MESH_TRANSFORM_PLATEMESH) { in NS_BLASIUS()
318 …PetscCall(ModifyMesh(comm, dm, problem->dim, mesh_growth, mesh_Ndelta, mesh_refine_height, mesh_to… in NS_BLASIUS()
323 …PetscCallCeed(ceed, CeedQFunctionContextGetData(problem->apply_vol_rhs.qfunction_context, CEED_MEM… in NS_BLASIUS()
331 blasius_ctx->weakT = weakT; in NS_BLASIUS()
332 blasius_ctx->T_wall = T_wall; in NS_BLASIUS()
333 blasius_ctx->delta0 = delta0; in NS_BLASIUS()
334 blasius_ctx->S_infty = S_infty; in NS_BLASIUS()
335 blasius_ctx->n_cheb = N; in NS_BLASIUS()
336 blasius_ctx->implicit = user->phys->implicit; in NS_BLASIUS()
337 …if (P0_set) newtonian_ig_ctx->idl_pressure = P_inf; // For maintaining behavior of -P0 flag (whic… in NS_BLASIUS()
338 blasius_ctx->newtonian_ctx = *newtonian_ig_ctx; in NS_BLASIUS()
343 blasius_ctx->x_inflow = domain_min[0]; in NS_BLASIUS()
344 blasius_ctx->eta_max = 5 * domain_max[1] / blasius_ctx->delta0; in NS_BLASIUS()
347 PetscCall(PetscOptionsGetBool(NULL, NULL, "-diff_filter_mms", &diff_filter_mms, NULL)); in NS_BLASIUS()
350 …PetscCallCeed(ceed, CeedQFunctionContextRestoreData(problem->apply_vol_rhs.qfunction_context, &new… in NS_BLASIUS()
352 PetscCallCeed(ceed, CeedQFunctionContextCreate(user->ceed, &blasius_context)); in NS_BLASIUS()
356 PetscCallCeed(ceed, CeedQFunctionContextDestroy(&problem->ics.qfunction_context)); in NS_BLASIUS()
357 problem->ics.qfunction_context = blasius_context; in NS_BLASIUS()
363 …PetscCheck((user->phys->state_var == STATEVAR_CONSERVATIVE) || (user->app_ctx->test_type == TESTTY… in NS_BLASIUS()
365 problem->apply_inflow.qfunction = Blasius_Inflow; in NS_BLASIUS()
366 problem->apply_inflow.qfunction_loc = Blasius_Inflow_loc; in NS_BLASIUS()
367 problem->apply_inflow_jacobian.qfunction = Blasius_Inflow_Jacobian; in NS_BLASIUS()
368 problem->apply_inflow_jacobian.qfunction_loc = Blasius_Inflow_Jacobian_loc; in NS_BLASIUS()
369 …PetscCallCeed(ceed, CeedQFunctionContextReferenceCopy(blasius_context, &problem->apply_inflow.qfun… in NS_BLASIUS()
370 …PetscCallCeed(ceed, CeedQFunctionContextReferenceCopy(blasius_context, &problem->apply_inflow_jaco… in NS_BLASIUS()