CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > General Forums > Main CFD Forum

On the Roe approximate Riemann solver in the preconditioned density based method

Register Blogs Community New Posts Updated Threads Search

Like Tree13Likes

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   March 14, 2022, 14:56
Default
  #21
Senior Member
 
Eifoehn4's Avatar
 
-
Join Date: Jul 2012
Location: Germany
Posts: 184
Rep Power: 14
Eifoehn4 is on a distinguished road
I don't know if this will help you:


Initially consider only the conservative Jacobian without preconditioner with

\boldsymbol{Q}=
   \begin{pmatrix}
      \rho      \\[0.2em]
      \rho u    \\[0.2em]
      \rho e  
  \end{pmatrix}=
   \begin{pmatrix}
      \alpha    \\[0.2em]
      \beta     \\[0.2em]
      \gamma  
  \end{pmatrix},
\quad \boldsymbol{F}(\boldsymbol{Q})=
   \begin{pmatrix}
      \rho u                                \\[0.2em]
      \rho u^2 + p                          \\[0.2em]
      u \left( E + p \right)    
  \end{pmatrix}=
   \begin{pmatrix}
     \beta                                  \\[0.2em]
      \beta^2/\alpha + p                    \\[0.2em]
      \beta \gamma/\alpha + p \beta/\alpha  
  \end{pmatrix},

and the pressure p=f(\alpha,\beta,\gamma) defined as a function of all conservative variables.

For non-ideal gases the crucial part is derive consistent Roe partial derivatives to fulfill following relation

\text{d}p=\left. \frac{\partial{p}}{\partial{\alpha}} \right|_{\beta,\gamma} \text{d}\alpha + \left. \frac{\partial{p}}{\partial{\beta}} \right|_{\alpha,\gamma} \text{d}\beta + \left. \frac{\partial{p}}{\partial{\gamma}} \right|_{\alpha,\beta} \text{d}\gamma.

One approach may be the Liou-Glaister generalization.

Here you initially calculate the partial derivatives point-wise with the EOS with the classical Roe averages (denoted with a tilde symbol).

They do not fulfill the first relation and result in a wrong \Delta p (origin of instability).

\Delta p \neq \tilde{p}_{\alpha}\Delta\alpha + \tilde{p}_{\beta}\Delta\beta + \tilde{p}_{\gamma}\Delta\gamma,

After that you correct them by a projection method (denoted with a hat symbol)

\partial p=- \Delta p + \tilde{p}_{\alpha}\Delta\alpha + \tilde{p}_{\beta}\Delta\beta + \tilde{p}_{\gamma}\Delta\gamma,

where

\hat{p}_{\alpha}    = \tilde{p}_{\alpha}    \left(1 + \tilde{p}_{\alpha}     \Delta \alpha \partial p / D \right) ,

\hat{p}_{\beta} = ...

\hat{p}_{\gamma}   = ...

with

D  = (\tilde{p}_{\alpha}\Delta\alpha)^2+(\tilde{p}_{\beta}\Delta\beta)^2+(\tilde{p}_{\gamma}\Delta\gamma)^2.

Now in analogy to that you may approach in a similar way for the preconditioner version by considering the two equations

\text{d}\rho=\left. \frac{\partial{\rho}}{\partial{p}} \right|_{T} \text{d}p + \left. \frac{\partial{\rho}}{\partial{T}} \right|_{p} \text{d}T + ...

and

\text{d}H=\left. \frac{\partial{H}}{\partial{p}} \right|_{T} \text{d}p + \left. \frac{\partial{H}}{\partial{T}} \right|_{p} \text{d}T + ...

to derive consistent primitive partial derivatives.

The approach is well described in the paper of Mottura.

Regards
sbaffini likes this.
__________________
Check out my side project:

A multiphysics discontinuous Galerkin framework: Youtube, Gitlab.

Last edited by Eifoehn4; March 15, 2022 at 14:59.
Eifoehn4 is offline   Reply With Quote

Old   March 15, 2022, 08:11
Default
  #22
Senior Member
 
sbaffini's Avatar
 
Paolo Lampitella
Join Date: Mar 2009
Location: Italy
Posts: 2,195
Blog Entries: 29
Rep Power: 39
sbaffini will become famous soon enoughsbaffini will become famous soon enough
Send a message via Skype™ to sbaffini
Quote:
Originally Posted by Eifoehn4 View Post
I don't know if this will help you:


Initially consider only the conservative Jacobian without preconditioner with

\boldsymbol{Q}=
   \begin{pmatrix}
      \rho      \\[0.2em]
      \rho u    \\[0.2em]
      \rho e  
  \end{pmatrix}=
   \begin{pmatrix}
      \alpha    \\[0.2em]
      \beta     \\[0.2em]
      \gamma  
  \end{pmatrix},
\quad \boldsymbol{F}(\boldsymbol{Q})=
   \begin{pmatrix}
      \rho u                                \\[0.2em]
      \rho u^2 + p                          \\[0.2em]
      u \left( E + p \right)    
  \end{pmatrix}=
   \begin{pmatrix}
     \beta                                  \\[0.2em]
      \beta^2/\alpha + p                    \\[0.2em]
      \beta \gamma/\alpha + p \beta/\alpha  
  \end{pmatrix},

and the pressure p=f(\alpha,\beta,\gamma) defined as a function of all conservative variables.

For non-ideal gases the crucial part is derive consistent Roe partial derivatives to fulfill following relation

\text{d}p=\left. \frac{\partial{p}}{\partial{\alpha}} \right|_{\beta,\gamma} \text{d}\alpha + \left. \frac{\partial{p}}{\partial{\beta}} \right|_{\alpha,\gamma} \text{d}\beta + \left. \frac{\partial{p}}{\partial{\gamma}} \right|_{\alpha,\beta} \text{d}\gamma.

One approach may be the Liou-Glaister generalization.

Here you initially calculate the partial derivatives point-wise with the EOS with the classical Roe averages (denoted with a tilde symbol).

\Delta p = \tilde{p}_{\alpha}\Delta\alpha + \tilde{p}_{\beta}\Delta\beta + \tilde{p}_{\gamma}\Delta\gamma,

They do not fulfill the first relation and result in a wrong \Delta p (origin of instability). After that you correct them by a projection method (denoted with a hat symbol)

\hat{p}_{\alpha}    = \tilde{p}_{\alpha}    \left(1 + \tilde{p}_{\alpha}     \Delta \alpha \Delta p / D \right) ,

\hat{p}_{\beta} = ...

\hat{p}_{\gamma}   = ...

with

D  = (\tilde{p}_{\alpha}\Delta\alpha)^2+(\tilde{p}_{\beta}\Delta\beta)^2+(\tilde{p}_{\gamma}\Delta\gamma)^2.

Now in analogy to that you may approach in a similar way for the preconditioner version by considering the two equations

\text{d}\rho=\left. \frac{\partial{\rho}}{\partial{p}} \right|_{T} \text{d}p + \left. \frac{\partial{\rho}}{\partial{T}} \right|_{p} \text{d}T + ...

and

\text{d}H=\left. \frac{\partial{H}}{\partial{p}} \right|_{T} \text{d}p + \left. \frac{\partial{H}}{\partial{T}} \right|_{p} \text{d}T + ...

to derive consistent primitive partial derivatives.

The approach is well described in the paper of Mottura.

Regards
I think I finally get how this should work. Still, I have a doubt about the fact that I would be applying a method intended to work for cases where the \delta X^* at numerator over D is necessarily not 0, to a case where, instead, it should be (but isn't, maybe only because I don't know how to write it with my variables). However, I might still use the idea to actually check what I'm doing.

For now, my euristic explanation for the matter is that, looking at my final dissipation vector here, my d1 component, which appears on all the equations, starts with drho/dT and all the other terms have rho_ROE inside, so drho/dT should sort of have it as well.
sbaffini is offline   Reply With Quote

Old   March 19, 2022, 12:20
Default
  #23
Senior Member
 
sbaffini's Avatar
 
Paolo Lampitella
Join Date: Mar 2009
Location: Italy
Posts: 2,195
Blog Entries: 29
Rep Power: 39
sbaffini will become famous soon enoughsbaffini will become famous soon enough
Send a message via Skype™ to sbaffini
Quote:
Originally Posted by sbaffini View Post
Just for your interest, this is what I have so far in my scheme for the dissipation part \mathbf{d}=\mathbf{\Gamma}\mathbf{R}|\mathbf{\Lambda}|\mathbf{R}^{-1}(\mathbf{q}_R-\mathbf{q}_L)(but note that I don't use it in this form, I still use matrix products to avoid loosing control):

d_1 = \rho_T|V_n|a+\rho b
d_2 = u d_1 + \rho|V_n|\left(u_R-u_L\right) + \rho n_x c
d_3 = v d_1 + \rho|V_n|\left(v_R-v_L\right) + \rho n_y c
d_4 = w d_1 + \rho|V_n|\left(w_R-w_L\right) + \rho n_z c
d_5 = H d_1 + \rho|V_n|e+ \rho V_n c

with:

a = \left(T_R-T_L\right)-\frac{\delta}{\rho C_p} \left(P_R-P_L\right)
b = \frac{\theta_2}{2\rho s_2V_r^2}\left(P_R-P_L\right)+\frac{\theta_3}{2s_2V_r^2}\left[\mathbf{n} \cdot \left(\mathbf{V}_R-\mathbf{V}_L\right)\right]
c = \frac{\theta_1}{2 \rho s_2} \left(P_R-P_L\right) +\left(\frac{\theta_2}{2 s_2} -|V_n| \right) \left[\mathbf{n} \cdot \left(\mathbf{V}_R-\mathbf{V}_L\right)\right]
e = C_p a + \mathbf{V} \cdot \left(\mathbf{V}_R-\mathbf{V}_L\right)

and:
\theta_0 = |s_1+s_2|+|s_1-s_2| = 2MAX\left(s1,s2\right)
\theta_1 = |s_1+s_2|-|s_1-s_2| = 2MIN\left(s1,s2\right)
\frac{\theta_2}{2s2} = \frac{\theta_1}{2s2}\left(V_n-s_1\right)-\frac{\theta_0}{2}
\frac{\theta_3}{2s2} = \frac{\theta_1}{2s2}\left[\left(V_n-s_1\right)^2+s_2^2\right]-\left(V_n-s_1\right)\theta_0
s_1 = \left(1-\alpha\right)V_n
s_2 = \sqrt{\left(\alpha V_n \right)^2+V_r^2}
\alpha = \frac{1-\beta V_r^2}{2}

where \beta=\rho_P+\rho T \delta /(\rho Cp) is the inverse of the speed of sound squared, \delta = 1 - \rho H_p=-\rho_T T/\rho and anything without the R/L subscript to be evaluated at the Roe average. One could actually further simplify things by nothing that:

V_n-s_1 = \alpha V_n
\left(V_n-s_1\right)^2+s_2^2=2 \alpha^2 V_n^2 + V_r^2

but I've found difficulties in proceeding further toward a version which could even remotely resemble the Roe one (not sure if it is also due to the certainly different eigenvectors I'm using).
Just a note to mention that the dissipation vector I reported above was actually wrong (as I inverted the 4th and 5th eigenvectors in doing the products). However, the error just affects the coefficients appearing in b and c, basically just \theta_2 and \theta_3. So, in the end, the heuristic reasoning on \rho_T still looks legit.
sbaffini is offline   Reply With Quote

Reply


Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
Pressure based and Density based Solver Xobile Main CFD Forum 39 August 19, 2020 07:04
Approximate Riemann solver and reconstruction TurbJet Main CFD Forum 2 June 14, 2020 10:02
Density and Pressure based solver Diger ANSYS 1 May 2, 2018 00:51
Density or Pressure based solver Achu FLUENT 1 April 27, 2018 05:34
Pressure based and Density based Solver taekyu8 FLUENT 0 January 28, 2013 12:05


All times are GMT -4. The time now is 00:05.