densitycurrent.c (f163c87f5f916c922650e26a0ad2b739121e1f13) densitycurrent.c (1864f1c2b4e770a2a9adc26a02ef77fc3a284256)
1// Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at
2// the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights
3// reserved. See files LICENSE and NOTICE for details.
4//
5// This file is part of CEED, a collection of benchmarks, miniapps, software
6// libraries and APIs for efficient high-order finite element and spectral
7// element discretizations for exascale applications. For more information and
8// source code availability see http://github.com/ceed.

--- 7 unchanged lines hidden (view full) ---

16
17/// @file
18/// Utility functions for setting up DENSITY_CURRENT
19
20#include "../navierstokes.h"
21#include "../qfunctions/setupgeo.h"
22#include "../qfunctions/densitycurrent.h"
23
1// Copyright (c) 2017, Lawrence Livermore National Security, LLC. Produced at
2// the Lawrence Livermore National Laboratory. LLNL-CODE-734707. All Rights
3// reserved. See files LICENSE and NOTICE for details.
4//
5// This file is part of CEED, a collection of benchmarks, miniapps, software
6// libraries and APIs for efficient high-order finite element and spectral
7// element discretizations for exascale applications. For more information and
8// source code availability see http://github.com/ceed.

--- 7 unchanged lines hidden (view full) ---

16
17/// @file
18/// Utility functions for setting up DENSITY_CURRENT
19
20#include "../navierstokes.h"
21#include "../qfunctions/setupgeo.h"
22#include "../qfunctions/densitycurrent.h"
23
24PetscErrorCode NS_DENSITY_CURRENT(ProblemData *problem, void *setup_ctx,
24PetscErrorCode NS_DENSITY_CURRENT(ProblemData *problem, DM dm, void *setup_ctx,
25 void *ctx) {
26 SetupContext setup_context = *(SetupContext *)setup_ctx;
27 User user = *(User *)ctx;
28 StabilizationType stab;
29 MPI_Comm comm = PETSC_COMM_WORLD;
30 PetscBool implicit;
31 PetscBool has_curr_time = PETSC_FALSE;
32 PetscInt ierr;

--- 34 unchanged lines hidden (view full) ---

67 CeedScalar cp = 1004.; // J/(kg K)
68 CeedScalar g = 9.81; // m/s^2
69 CeedScalar lambda = -2./3.; // -
70 CeedScalar mu = 75.; // Pa s, dynamic viscosity
71 // mu = 75 is not physical for air, but is good for numerical stability
72 CeedScalar k = 0.02638; // W/(m K)
73 CeedScalar c_tau = 0.5; // -
74 // c_tau = 0.5 is reported as "optimal" in Hughes et al 2010
25 void *ctx) {
26 SetupContext setup_context = *(SetupContext *)setup_ctx;
27 User user = *(User *)ctx;
28 StabilizationType stab;
29 MPI_Comm comm = PETSC_COMM_WORLD;
30 PetscBool implicit;
31 PetscBool has_curr_time = PETSC_FALSE;
32 PetscInt ierr;

--- 34 unchanged lines hidden (view full) ---

67 CeedScalar cp = 1004.; // J/(kg K)
68 CeedScalar g = 9.81; // m/s^2
69 CeedScalar lambda = -2./3.; // -
70 CeedScalar mu = 75.; // Pa s, dynamic viscosity
71 // mu = 75 is not physical for air, but is good for numerical stability
72 CeedScalar k = 0.02638; // W/(m K)
73 CeedScalar c_tau = 0.5; // -
74 // c_tau = 0.5 is reported as "optimal" in Hughes et al 2010
75 PetscScalar lx = 8000.; // m
76 PetscScalar ly = 8000.; // m
77 PetscScalar lz = 4000.; // m
78 CeedScalar rc = 1000.; // m (Radius of bubble)
79 PetscReal center[3], dc_axis[3] = {0, 0, 0};
75 CeedScalar rc = 1000.; // m (Radius of bubble)
76 PetscReal center[3], dc_axis[3] = {0, 0, 0};
77 PetscReal domain_min[3], domain_max[3], domain_size[3];
78 ierr = DMGetBoundingBox(dm, domain_min, domain_max); CHKERRQ(ierr);
79 for (int i=0; i<3; i++) domain_size[i] = domain_max[i] - domain_min[i];
80
81 // ------------------------------------------------------
82 // Create the PETSc context
83 // ------------------------------------------------------
84 PetscScalar meter = 1e-2; // 1 meter in scaled length units
85 PetscScalar kilogram = 1e-6; // 1 kilogram in scaled mass units
86 PetscScalar second = 1e-2; // 1 second in scaled time units
87 PetscScalar Kelvin = 1; // 1 Kelvin in scaled temperature units

--- 21 unchanged lines hidden (view full) ---

109 NULL, g, &g, NULL); CHKERRQ(ierr);
110 ierr = PetscOptionsScalar("-lambda",
111 "Stokes hypothesis second viscosity coefficient",
112 NULL, lambda, &lambda, NULL); CHKERRQ(ierr);
113 ierr = PetscOptionsScalar("-mu", "Shear dynamic viscosity coefficient",
114 NULL, mu, &mu, NULL); CHKERRQ(ierr);
115 ierr = PetscOptionsScalar("-k", "Thermal conductivity",
116 NULL, k, &k, NULL); CHKERRQ(ierr);
80
81 // ------------------------------------------------------
82 // Create the PETSc context
83 // ------------------------------------------------------
84 PetscScalar meter = 1e-2; // 1 meter in scaled length units
85 PetscScalar kilogram = 1e-6; // 1 kilogram in scaled mass units
86 PetscScalar second = 1e-2; // 1 second in scaled time units
87 PetscScalar Kelvin = 1; // 1 Kelvin in scaled temperature units

--- 21 unchanged lines hidden (view full) ---

109 NULL, g, &g, NULL); CHKERRQ(ierr);
110 ierr = PetscOptionsScalar("-lambda",
111 "Stokes hypothesis second viscosity coefficient",
112 NULL, lambda, &lambda, NULL); CHKERRQ(ierr);
113 ierr = PetscOptionsScalar("-mu", "Shear dynamic viscosity coefficient",
114 NULL, mu, &mu, NULL); CHKERRQ(ierr);
115 ierr = PetscOptionsScalar("-k", "Thermal conductivity",
116 NULL, k, &k, NULL); CHKERRQ(ierr);
117 ierr = PetscOptionsScalar("-lx", "Length scale in x direction",
118 NULL, lx, &lx, NULL); CHKERRQ(ierr);
119 ierr = PetscOptionsScalar("-ly", "Length scale in y direction",
120 NULL, ly, &ly, NULL); CHKERRQ(ierr);
121 ierr = PetscOptionsScalar("-lz", "Length scale in z direction",
122 NULL, lz, &lz, NULL); CHKERRQ(ierr);
123 ierr = PetscOptionsScalar("-rc", "Characteristic radius of thermal bubble",
124 NULL, rc, &rc, NULL); CHKERRQ(ierr);
117 ierr = PetscOptionsScalar("-rc", "Characteristic radius of thermal bubble",
118 NULL, rc, &rc, NULL); CHKERRQ(ierr);
119 for (int i=0; i<3; i++) center[i] = .5*domain_size[i];
125 PetscInt n = problem->dim;
120 PetscInt n = problem->dim;
126 center[0] = 0.5 * lx;
127 center[1] = 0.5 * ly;
128 center[2] = 0.5 * lz;
129 ierr = PetscOptionsRealArray("-center", "Location of bubble center",
130 NULL, center, &n, NULL); CHKERRQ(ierr);
131 n = problem->dim;
132 ierr = PetscOptionsRealArray("-dc_axis",
133 "Axis of density current cylindrical anomaly, or {0,0,0} for spherically symmetric",
134 NULL, dc_axis, &n, NULL); CHKERRQ(ierr);
135 {
136 PetscReal norm = PetscSqrtReal(PetscSqr(dc_axis[0]) + PetscSqr(dc_axis[1]) +

--- 61 unchanged lines hidden (view full) ---

198 thetaC *= Kelvin;
199 P0 *= Pascal;
200 N *= (1./second);
201 cv *= J_per_kg_K;
202 cp *= J_per_kg_K;
203 g *= m_per_squared_s;
204 mu *= Pascal * second;
205 k *= W_per_m_K;
121 ierr = PetscOptionsRealArray("-center", "Location of bubble center",
122 NULL, center, &n, NULL); CHKERRQ(ierr);
123 n = problem->dim;
124 ierr = PetscOptionsRealArray("-dc_axis",
125 "Axis of density current cylindrical anomaly, or {0,0,0} for spherically symmetric",
126 NULL, dc_axis, &n, NULL); CHKERRQ(ierr);
127 {
128 PetscReal norm = PetscSqrtReal(PetscSqr(dc_axis[0]) + PetscSqr(dc_axis[1]) +

--- 61 unchanged lines hidden (view full) ---

190 thetaC *= Kelvin;
191 P0 *= Pascal;
192 N *= (1./second);
193 cv *= J_per_kg_K;
194 cp *= J_per_kg_K;
195 g *= m_per_squared_s;
196 mu *= Pascal * second;
197 k *= W_per_m_K;
206 lx = fabs(lx) * meter;
207 ly = fabs(ly) * meter;
208 lz = fabs(lz) * meter;
209 rc = fabs(rc) * meter;
198 rc = fabs(rc) * meter;
199 for (int i=0; i<3; i++) domain_size[i] *= meter;
210 for (int i=0; i<3; i++) center[i] *= meter;
200 for (int i=0; i<3; i++) center[i] *= meter;
201 problem->dm_scale = meter;
211
212 // -- Setup Context
213 setup_context->theta0 = theta0;
214 setup_context->thetaC = thetaC;
215 setup_context->P0 = P0;
216 setup_context->N = N;
217 setup_context->cv = cv;
218 setup_context->cp = cp;
219 setup_context->g = g;
220 setup_context->rc = rc;
202
203 // -- Setup Context
204 setup_context->theta0 = theta0;
205 setup_context->thetaC = thetaC;
206 setup_context->P0 = P0;
207 setup_context->N = N;
208 setup_context->cv = cv;
209 setup_context->cp = cp;
210 setup_context->g = g;
211 setup_context->rc = rc;
221 setup_context->lx = lx;
222 setup_context->ly = ly;
223 setup_context->lz = lz;
212 setup_context->lx = domain_size[0];
213 setup_context->ly = domain_size[1];
214 setup_context->lz = domain_size[2];
224 setup_context->center[0] = center[0];
225 setup_context->center[1] = center[1];
226 setup_context->center[2] = center[2];
227 setup_context->dc_axis[0] = dc_axis[0];
228 setup_context->dc_axis[1] = dc_axis[1];
229 setup_context->dc_axis[2] = dc_axis[2];
230 setup_context->time = 0;
231

--- 147 unchanged lines hidden --- 		215 setup_context->center[0] = center[0];
216 setup_context->center[1] = center[1];
217 setup_context->center[2] = center[2];
218 setup_context->dc_axis[0] = dc_axis[0];
219 setup_context->dc_axis[1] = dc_axis[1];
220 setup_context->dc_axis[2] = dc_axis[2];
221 setup_context->time = 0;
222

--- 147 unchanged lines hidden ---