1 // Copyright (c) 2017-2022, Lawrence Livermore National Security, LLC and other CEED contributors. 2 // All Rights Reserved. See the top-level LICENSE and NOTICE files for details. 3 // 4 // SPDX-License-Identifier: BSD-2-Clause 5 // 6 // This file is part of CEED: http://github.com/ceed 7 8 /// @file 9 /// Utility functions for setting up DENSITY_CURRENT 10 11 #include "../navierstokes.h" 12 #include "../qfunctions/densitycurrent.h" 13 14 PetscErrorCode NS_DENSITY_CURRENT(ProblemData *problem, DM dm, void *setup_ctx, 15 void *ctx) { 16 17 PetscInt ierr; 18 ierr = NS_NEWTONIAN_IG(problem, dm, setup_ctx, ctx); CHKERRQ(ierr); 19 SetupContext setup_context = *(SetupContext *)setup_ctx; 20 User user = *(User *)ctx; 21 MPI_Comm comm = PETSC_COMM_WORLD; 22 PetscFunctionBeginUser; 23 24 // ------------------------------------------------------ 25 // SET UP DENSITY_CURRENT 26 // ------------------------------------------------------ 27 problem->ics = ICsDC; 28 problem->ics_loc = ICsDC_loc; 29 problem->bc = Exact_DC; 30 31 // ------------------------------------------------------ 32 // Create the libCEED context 33 // ------------------------------------------------------ 34 CeedScalar rc = 1000.; // m (Radius of bubble) 35 PetscReal center[3], dc_axis[3] = {0, 0, 0}; 36 PetscReal domain_min[3], domain_max[3], domain_size[3]; 37 ierr = DMGetBoundingBox(dm, domain_min, domain_max); CHKERRQ(ierr); 38 for (int i=0; i<3; i++) domain_size[i] = domain_max[i] - domain_min[i]; 39 40 // ------------------------------------------------------ 41 // Command line Options 42 // ------------------------------------------------------ 43 ierr = PetscOptionsBegin(comm, NULL, "Options for DENSITY_CURRENT problem", 44 NULL); CHKERRQ(ierr); 45 ierr = PetscOptionsScalar("-rc", "Characteristic radius of thermal bubble", 46 NULL, rc, &rc, NULL); CHKERRQ(ierr); 47 for (int i=0; i<3; i++) center[i] = .5*domain_size[i]; 48 PetscInt n = problem->dim; 49 ierr = PetscOptionsRealArray("-center", "Location of bubble center", 50 NULL, center, &n, NULL); CHKERRQ(ierr); 51 n = problem->dim; 52 ierr = PetscOptionsRealArray("-dc_axis", 53 "Axis of density current cylindrical anomaly, or {0,0,0} for spherically symmetric", 54 NULL, dc_axis, &n, NULL); CHKERRQ(ierr); 55 { 56 PetscReal norm = PetscSqrtReal(PetscSqr(dc_axis[0]) + PetscSqr(dc_axis[1]) + 57 PetscSqr(dc_axis[2])); 58 if (norm > 0) { 59 for (int i=0; i<3; i++) dc_axis[i] /= norm; 60 } 61 } 62 63 ierr = PetscOptionsEnd(); CHKERRQ(ierr); 64 65 PetscScalar meter = user->units->meter; 66 rc = fabs(rc) * meter; 67 for (int i=0; i<3; i++) center[i] *= meter; 68 69 setup_context->rc = rc; 70 setup_context->center[0] = center[0]; 71 setup_context->center[1] = center[1]; 72 setup_context->center[2] = center[2]; 73 setup_context->dc_axis[0] = dc_axis[0]; 74 setup_context->dc_axis[1] = dc_axis[1]; 75 setup_context->dc_axis[2] = dc_axis[2]; 76 77 PetscFunctionReturn(0); 78 } 79 80 PetscErrorCode SetupContext_DENSITY_CURRENT(Ceed ceed, CeedData ceed_data, 81 AppCtx app_ctx, SetupContext setup_ctx, Physics phys) { 82 PetscFunctionBeginUser; 83 PetscInt ierr = SetupContext_NEWTONIAN_IG(ceed, ceed_data, app_ctx, setup_ctx, 84 phys); 85 CHKERRQ(ierr); 86 PetscFunctionReturn(0); 87 } 88 89 PetscErrorCode PRINT_DENSITY_CURRENT(Physics phys, SetupContext setup_ctx, 90 AppCtx app_ctx) { 91 MPI_Comm comm = PETSC_COMM_WORLD; 92 PetscErrorCode ierr; 93 PetscFunctionBeginUser; 94 95 ierr = PetscPrintf(comm, 96 " Problem:\n" 97 " Problem Name : %s\n" 98 " Stabilization : %s\n", 99 app_ctx->problem_name, StabilizationTypes[phys->stab]); 100 CHKERRQ(ierr); 101 102 PetscFunctionReturn(0); 103 } 104