I am trying to implement a 1D k-$\varepsilon$ turbulent model for a turbidity current, hence the conservation equation for $c$. I'm solving for the variables $u,c,k$ and $\varepsilon$. The remaining terms are parameters. $$ \newcommand{\cmu}{c_{\mu}} \newcommand{\cmuc}{c_{\mu c}} \newcommand{\dervz}[1]{\frac{\partial #1}{\partial z}} \newcommand{\eddyvis}{\cmu \frac{k^2}{\varepsilon}} \newcommand{\eddyvisC}{\cmuc \frac{\hat{k}^2}{\hat{\varepsilon}}} $$

$$ \frac{\partial u }{\partial t} = \frac{\partial}{\partial z} \left[ \left( \nu + \cmu \frac{k^2}{\varepsilon} \right) \frac{\partial u }{\partial z} \right] + \frac{\rho_s - \rho_w}{\rho_w} c g_x $$

$$ \frac{\partial c}{\partial t} = v_s \dervz{c} + \frac{\partial}{\partial z} \left[ \left( \nu + \cmuc \frac{k^2}{\varepsilon} \right) \dervz{c} \right] $$

$$ \frac{\partial k}{\partial t} = \frac{\partial}{\partial z} \left[ \left( \nu + \frac{\cmu}{\sigma_k} \frac{k^2}{\varepsilon} \right) \dervz{k} \right] + \eddyvis \left(\dervz{U}\right)^2 - \varepsilon + \frac{\rho - \rho_w}{\rho_w} g_z \frac{1}{\sigma} \eddyvisC \dervz{c} $$

$$ \frac{\partial \varepsilon}{\partial t} = \frac{\partial}{\partial z} \left[ \left( \nu + \frac{\cmu}{\sigma_{\varepsilon}} \frac{k^2}{\varepsilon} \right) \dervz{\varepsilon} \right] + c_{\varepsilon 1} \frac{\varepsilon}{k} \left( \eddyvis \left(\dervz{u}\right)^2 + c_{\varepsilon 3} \frac{\rho - \rho_w}{\rho_w} g_z \frac{1}{\sigma} \eddyvisC \dervz{c} \right) - C_{\varepsilon 2} \frac{ \varepsilon^2}{k} $$

So far, I am using finite difference. Upwind scheme for the first order derivatives. Because the paramter $v_s$ is positive, concentration travels downwards at a speed $v_s$, I'm implementing the upwind scheme as: $$ \frac{c^{i+1}_n - c^{i}_n}{\Delta z} $$ for the node $i$. I also use the same upwind scheme for all the first-order derivative terms, i.e. also the terms inside the diffusion terms, which are nonlinear. For the linear diffusion terms, I apply second order central difference. For the time integration, I use a Crank-Nicholson method treating all terms implicitly.

The boundary conditions are complex. I'm using a standard wall function, using the log law at a node far from the wall. Using the log law, I also impose a gradient boundary condition because the log law adds the shear velocity as a unknown, hence the additional equation/condition. For $k$ and $\varepsilon$ I use conditions that depend on the shear velocity. That's for one side of the domain, at the other side, sufficiently far, I use zero value conditions.

I'm not seeing much success with this implementation, should I use different discretization schemes for the diffusion and the advective terms? Should I use a different method altogether such as Finite Volume? What confuses me are the diffusion terms. Which discretization should I use there?

With regards to the time discretization. Should I treat some terms explicitely? The source terms are nonlinear as well.

Thanks in advance.

EDIT: Response to the first comment

Sorry I was not clear there. The model blows up for certain initial values of the concentration (too much acceleration), although that might have to do with the wall function implementation I think. Another issue that I had were the initial conditions for $k$ and $\varepsilon$. I'll explain. The initial conditions for the velocity and concentration do not span the entire domain. They have positive values only for one tenth of the lower part of the domain, where I implement the wall function boundary condition, eventually, they will diffuse upwards. However, they will not do so if the values for $k$ and $\varepsilon$ are zero above where $u$ and $c$ are. This is because they are in the diffusion term and I think because of the upwind scheme, which assumes the wave travels downwards.

  • $\begingroup$ What are you out for? Since turbulence is a 3D phenomenon, I don't think that $\kappa-\epsilon$ models make any sense in 1D. $\endgroup$ – Jan Jul 19 '15 at 18:18
  • $\begingroup$ It's a one dimensional velocity profile, where we have neglected the velocity components that are not the streamwise components. There's still turbulence since there are velocities gradients. $\endgroup$ – balborian Jul 19 '15 at 19:38
  • $\begingroup$ I don't think this is enough for turbulence (models) to have a velocity gradient only in streamwise direction. E.g., the basic models for turbulence, like the $\kappa-\epsilon$ model, base on Prandtl's mixing length model, that bases on $\partial \bar u /\partial y$, i.e., the gradient of the mean flow in transversal direction. $\endgroup$ – Jan Jul 20 '15 at 6:34
  • $\begingroup$ There are no velocity component in the transversal direction, but there are velocity gradients. I typed them as $\frac{\partial u}{\partial z}$. $\endgroup$ – balborian Jul 20 '15 at 13:59
  • $\begingroup$ USe crank nicolson coz its more stable and conserves mass $\endgroup$ – Royston Fernandes Nov 6 '17 at 10:06

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.