what's the meaning of UEqn().A()

wangle · May 2, 2010, 10:59

I'm a new user of OpenFOAM,now I'm analyzing a solver .But I don't know what's the meaning of UEqn().A() and UEqn().H().Can you help me translate them to mathematic form.Thank you!

CedricVH · May 2, 2010, 11:13

Terms dependent on U are included in the A field.
Terms not directly dependent on U are included in the H field.

wangle · May 2, 2010, 11:25

Quote:

Originally Posted by CedricVH

Terms dependent on U are included in the A field.
Terms not directly dependent on U are included in the H field.

Thank you for your reply.Does it mean the matrix A in equntion AU=H?But it's expressed a volScalarField，what's the connection between the matrix and the "scalar"? Thank you!

mchurchf · May 3, 2010, 11:21

In the solver you are analyzing, the equation system to be solved begins as CU = R where C is a matrix, U is the solution vector, and R is the right hand side. The C matrix can be split into a matrix with only the diagonal elements of C, which is called A, and a matrix that has only the off-diagonal elements of C, which is called H'. In other words C = A + H'.

Therefore, the linear system becomes (A + H')U = R, which is the same as AU = R - H'U. The right hand side is simply called H, so H = R - H'U. Therefore, we know have AU = H.

So take a look at the code you are analyzing, and you'll see something similar to:

Code:

    fvVectorMatrix UEqn
    (
        fvm::ddt(U)                        // time derivative
      + fvm::div(phi, U)                  // convection
      + turbulence->divDevReff(U)    // viscous and turbulent deviatoric stresses
      ==
      - gradPd                               // specified mean pressure gradient
    );

Later you'll see code like rUA = 1.0/UEqn.A() and U = rUA*UEqn.H(). The .A() operator gives A formed in UEqn, the list of diagonal elements as explained above. This is a scalarField because each element is a scalar and corresponds to one grid cell. The .H() operator is a list of the elements of the vector H described above and formed in UEqn. If the variable U is a vector, then each element of H will be a vector; if the variable U is a scalar, then each element of H will be a vector, and so on.

It is important to note, though, that there will be a piece of code that forms an fv<Type>Matrix, like the piece of code I included above. Later, there will be a piece of code that says something like

Code:

solve(UEqn == - fvc::grad(pd) - fvc::grad(rhok) * gh);

where the matrix system is actually solved. Note that it is UEqn == ... The terms on the right hand side of the == do not contribute the the A and H given by UEqn.A() and UEqn.H(). Only the terms contained inside UEqn definition itself create UEqn.A() and UEqn.H().

For more description of this go take at look at section 2 of my description of buoyantBoussinesqPisoFoam at http://openfoamwiki.net/index.php/Bu...sinesqPisoFoam

wangle · May 4, 2010, 04:36

Quote:

Originally Posted by mchurchf

In the solver you are analyzing, the equation system to be solved begins as CU = R where C is a matrix, U is the solution vector, and R is the right hand side. The C matrix can be split into a matrix with only the diagonal elements of C, which is called A, and a matrix that has only the off-diagonal elements of C, which is called H'. In other words C = A + H'.

Therefore, the linear system becomes (A + H')U = R, which is the same as AU = R - H'U. The right hand side is simply called H, so H = R - H'U. Therefore, we know have AU = H.

So take a look at the code you are analyzing, and you'll see something similar to:

Code:

    fvVectorMatrix UEqn
    (
        fvm::ddt(U)                        // time derivative
      + fvm::div(phi, U)                  // convection
      + turbulence->divDevReff(U)    // viscous and turbulent deviatoric stresses
      ==
      - gradPd                               // specified mean pressure gradient
    );

Later you'll see code like rUA = 1.0/UEqn.A() and U = rUA*UEqn.H(). The .A() operator gives A formed in UEqn, the list of diagonal elements as explained above. This is a scalarField because each element is a scalar and corresponds to one grid cell. The .H() operator is a list of the elements of the vector H described above and formed in UEqn. If the variable U is a vector, then each element of H will be a vector; if the variable U is a scalar, then each element of H will be a vector, and so on.

It is important to note, though, that there will be a piece of code that forms an fv<Type>Matrix, like the piece of code I included above. Later, there will be a piece of code that says something like

Code:

solve(UEqn == - fvc::grad(pd) - fvc::grad(rhok) * gh);

where the matrix system is actually solved. Note that it is UEqn == ... The terms on the right hand side of the == do not contribute the the A and H given by UEqn.A() and UEqn.H(). Only the terms contained inside UEqn definition itself create UEqn.A() and UEqn.H().

For more description of this go take at look at section 2 of my description of buoyantBoussinesqPisoFoam at http://openfoamwiki.net/index.php/Bu...sinesqPisoFoam

Thank you very much to help me solve the problem!

beauty · May 6, 2010, 03:23

Quote:

Originally Posted by mchurchf

For more description of this go take at look at section 2 of my description of buoyantBoussinesqPisoFoam at http://openfoamwiki.net/index.php/Bu...sinesqPisoFoam

Hi, mchurchf
After reading your description of buoyantBoussinesPisoFoam, I am still puzzled about the meaning of the red line in the following code, which is in the UEqn of twophaseEulerFoam.I can not find the corresponding mathematical expression, can you help me? Thank you!
UaEqn =
(
(scalar(1) + Cvm*rhob*beta/rhoa)*
(
fvm::ddt(Ua)
+ fvm::div(phia, Ua, "div(phia,Ua)")
- fvm::Sp(fvc::div(phia), Ua)
)
- fvm::laplacian(nuEffa, Ua)
+ fvc::div(Rca)
+ fvm::div(phiRa, Ua, "div(phia,Ua)")
- fvm::Sp(fvc::div(phiRa), Ua)
+ (fvc::grad(alpha)/(fvc::average(alpha) + scalar(0.001)) & Rca)
==
- fvm::Sp(beta/rhoa*K, Ua)
- beta/rhoa*(liftCoeff - Cvm*rhob*DDtUb)

Duo · November 9, 2016, 09:22

Quote:

Originally Posted by wangle

I'm a new user of OpenFOAM,now I'm analyzing a solver .But I don't know what's the meaning of UEqn().A() and UEqn().H().Can you help me translate them to mathematic form.Thank you!

I am a new user too, I wonder how to know there is a member function named A() or H() of UEqn? where can I find the definition of UEqn, I checked UEqn.H, but I can't find it, anyone can give a link of the definition of some instruction? thanks!

khedar · November 9, 2016, 11:50

Hi Duo,
UEqn is and object of class fvVectorMatrix so you need to search for member functions of this class to get A() or U().

Also fvVectorMatrix is basically typedef name for fvMatrix<vector>.

Code:

typedef fvMatrix<vector> fvVectorMatrix

typedef fvMatrix<vector> fvVectorMatrix

so you search for fvMatrix<Type> template class for finding the required member functions. This can be found here:
http://openfoam.com/documentation/cp...ml/a00955.html

Duo · November 9, 2016, 11:55

thanks, very helpful!

usv001 · April 27, 2018, 04:38

Hello FOAMers,

I am still not sure about one point in Matthew's explanation:

Quote:

Originally Posted by mchurchf

The .A() operator gives A formed in UEqn, the list of diagonal elements as explained above. This is a scalarField because each element is a scalar and corresponds to one grid cell.

This may seem like a silly question but since UEqn is a vector equation, there should be three matrices (one each for u, v and w) and, correspondingly, three diagonal coefficients. So, how come the '.A()' operator gives a scalarField instead of a vectorField?

Many thanks.

Regards,
USV

usv001 · May 2, 2018, 08:09

Dear FOAMers,

After a quick derivation, I found out that all coefficients of the UEqn should be scalars. I am providing the derivation here for other users who may have similar doubts.

Consider the discretization of the convective term $\nabla\cdot \left(\bold{uu}\right)$ on a non-uniform 2D Cartesian grid shown below:

Code:

    -------
    |  N  |
---------------
| W |  P  | E |
---------------
    |  S  |
    -------

Integrating $\nabla\cdot \left(\bold{uu}\right)$ over cell P and applying Gauss theorem, we get:

$\iint_{\Omega_P}{\nabla\cdot \left(\bold{uu}\right)}\ dV = \oint_{\partial \Omega_P}{\left(\bold{uu}\right)\cdot \mathbf{dS}} \approx \sum_{f}{\overbrace{\left(\mathbf{u}^0_f\cdot \mathbf{S}_f \right)}^{\phi^0_f} \mathbf{u}_f}$

If we use the central (linear) scheme to estimate $\phi^0_f$ and $\mathbf{u}_f$ at each face, we obtain:

Face e
$\phi_e^0 = \left[ \left(\frac{x_e - x_P}{x_E - x_P}\right)u_E^0 + \left(\frac{x_E - x_e}{x_E - x_P}\right)u_P^0 \right] \Delta{}y_e$
$\mathbf{u}_e = \left[ \left(\frac{x_e - x_P}{x_E - x_P}\right)\mathbf{u}_E + \left(\frac{x_E - x_e}{x_E - x_P}\right)\mathbf{u}_P \right]$

Face w
$\phi_w^0 = \left[ \left(\frac{x_w - x_P}{x_W - x_P}\right)u_W^0 + \left(\frac{x_W - x_w}{x_W - x_P}\right)u_P^0 \right] \Delta{}y_w$
$\mathbf{u}_w = \left[ \left(\frac{x_w - x_P}{x_W - x_P}\right)\mathbf{u}_W + \left(\frac{x_W - x_w}{x_W - x_P}\right)\mathbf{u}_P \right]$

Face n
$\phi_n^0 = \left[ \left(\frac{y_n - y_P}{y_N - y_P}\right)v_N^0 + \left(\frac{y_N - y_n}{y_N - y_P}\right)v_P^0 \right] \Delta{}x_n$
$\mathbf{u}_n = \left[ \left(\frac{y_n - y_P}{y_N - y_P}\right)\mathbf{u}_N + \left(\frac{y_N - y_n}{y_N - y_P}\right)\mathbf{u}_P \right]$

Face s
$\phi_s^0 = \left[ \left(\frac{y_s - y_P}{y_S - y_P}\right)v_S^0 + \left(\frac{y_S - y_s}{y_S - y_P}\right)v_P^0 \right] \Delta{}x_s$
$\mathbf{u}_s = \left[ \left(\frac{y_s - y_P}{y_S - y_P}\right)\mathbf{u}_S + \left(\frac{y_S - y_s}{y_S - y_P}\right)\mathbf{u}_P \right]$

Notice that all the coefficients are scalars. Now, these can be assembled in the form shown below. The rest of the discussion in this thread follows from here.

$a_P \mathbf{u}_P + \sum_{nb}{a_{nb} \mathbf{u}_{nb}}$

Hope this helps.

Cheers,
USV

May 2, 2010, 10:59	what's the meaning of UEqn().A()	#1
wangle New Member wangle Join Date: Dec 2009 Posts: 7 Rep Power: 17	I'm a new user of OpenFOAM,now I'm analyzing a solver .But I don't know what's the meaning of UEqn().A() and UEqn().H().Can you help me translate them to mathematic form.Thank you! Kummi and nepomnyi like this.

May 2, 2010, 11:13		#2
CedricVH Member Cedric Van Holsbeke Join Date: Dec 2009 Location: Belgium Posts: 81 Rep Power: 17	Terms dependent on U are included in the A field. Terms not directly dependent on U are included in the H field.

May 3, 2010, 11:21		#4
mchurchf Member Matthew J. Churchfield Join Date: Nov 2009 Location: Boulder, Colorado, USA Posts: 49 Rep Power: 19	In the solver you are analyzing, the equation system to be solved begins as CU = R where C is a matrix, U is the solution vector, and R is the right hand side. The C matrix can be split into a matrix with only the diagonal elements of C, which is called A, and a matrix that has only the off-diagonal elements of C, which is called H'. In other words C = A + H'. Therefore, the linear system becomes (A + H')U = R, which is the same as AU = R - H'U. The right hand side is simply called H, so H = R - H'U. Therefore, we know have AU = H. So take a look at the code you are analyzing, and you'll see something similar to: Code: fvVectorMatrix UEqn ( fvm::ddt(U) // time derivative + fvm::div(phi, U) // convection + turbulence->divDevReff(U) // viscous and turbulent deviatoric stresses == - gradPd // specified mean pressure gradient ); Later you'll see code like rUA = 1.0/UEqn.A() and U = rUAUEqn.H(). The .A() operator gives A formed in UEqn, the list of diagonal elements as explained above. This is a scalarField because each element is a scalar and corresponds to one grid cell. The .H() operator is a list of the elements of the vector H described above and formed in UEqn. If the variable U is a vector, then each element of H will be a vector; if the variable U is a scalar, then each element of H will be a vector, and so on. It is important to note, though, that there will be a piece of code that forms an fv<Type>Matrix, like the piece of code I included above. Later, there will be a piece of code that says something like Code: solve(UEqn == - fvc::grad(pd) - fvc::grad(rhok) gh); where the matrix system is actually solved. Note that it is UEqn == ... The terms on the right hand side of the == do not contribute the the A and H given by UEqn.A() and UEqn.H(). Only the terms contained inside UEqn definition itself create UEqn.A() and UEqn.H(). For more description of this go take at look at section 2 of my description of buoyantBoussinesqPisoFoam at http://openfoamwiki.net/index.php/Bu...sinesqPisoFoam karamiag, bigphil, Pascal_doran and 75 others like this.

November 9, 2016, 11:50		#8
khedar Senior Member khedar Join Date: Oct 2016 Posts: 111 Rep Power: 10	Hi Duo, UEqn is and object of class fvVectorMatrix so you need to search for member functions of this class to get A() or U(). Also fvVectorMatrix is basically typedef name for fvMatrix<vector>. Code: typedef fvMatrix<vector> fvVectorMatrix typedef fvMatrix<vector> fvVectorMatrix so you search for fvMatrix<Type> template class for finding the required member functions. This can be found here: http://openfoam.com/documentation/cp...ml/a00955.html acs, Duo and Kummi like this.

May 2, 2018, 08:09		#11
usv001 Senior Member Join Date: Sep 2015 Location: Singapore Posts: 102 Rep Power: 11	Dear FOAMers, After a quick derivation, I found out that all coefficients of the UEqn should be scalars. I am providing the derivation here for other users who may have similar doubts. Consider the discretization of the convective term $\nabla\cdot \left(\bold{uu}\right)$ on a non-uniform 2D Cartesian grid shown below: Code: ------- \| N \| --------------- \| W \| P \| E \| --------------- \| S \| ------- Integrating $\nabla\cdot \left(\bold{uu}\right)$ over cell P and applying Gauss theorem, we get: $\iint_{\Omega_P}{\nabla\cdot \left(\bold{uu}\right)}\ dV = \oint_{\partial \Omega_P}{\left(\bold{uu}\right)\cdot \mathbf{dS}} \approx \sum_{f}{\overbrace{\left(\mathbf{u}^0_f\cdot \mathbf{S}_f \right)}^{\phi^0_f} \mathbf{u}_f}$ If we use the central (linear) scheme to estimate $\phi^0_f$ and $\mathbf{u}_f$ at each face, we obtain: Face e $\phi_e^0 = \left[ \left(\frac{x_e - x_P}{x_E - x_P}\right)u_E^0 + \left(\frac{x_E - x_e}{x_E - x_P}\right)u_P^0 \right] \Delta{}y_e$ $\mathbf{u}_e = \left[ \left(\frac{x_e - x_P}{x_E - x_P}\right)\mathbf{u}_E + \left(\frac{x_E - x_e}{x_E - x_P}\right)\mathbf{u}_P \right]$ Face w $\phi_w^0 = \left[ \left(\frac{x_w - x_P}{x_W - x_P}\right)u_W^0 + \left(\frac{x_W - x_w}{x_W - x_P}\right)u_P^0 \right] \Delta{}y_w$ $\mathbf{u}_w = \left[ \left(\frac{x_w - x_P}{x_W - x_P}\right)\mathbf{u}_W + \left(\frac{x_W - x_w}{x_W - x_P}\right)\mathbf{u}_P \right]$ Face n $\phi_n^0 = \left[ \left(\frac{y_n - y_P}{y_N - y_P}\right)v_N^0 + \left(\frac{y_N - y_n}{y_N - y_P}\right)v_P^0 \right] \Delta{}x_n$ $\mathbf{u}_n = \left[ \left(\frac{y_n - y_P}{y_N - y_P}\right)\mathbf{u}_N + \left(\frac{y_N - y_n}{y_N - y_P}\right)\mathbf{u}_P \right]$ Face s $\phi_s^0 = \left[ \left(\frac{y_s - y_P}{y_S - y_P}\right)v_S^0 + \left(\frac{y_S - y_s}{y_S - y_P}\right)v_P^0 \right] \Delta{}x_s$ $\mathbf{u}_s = \left[ \left(\frac{y_s - y_P}{y_S - y_P}\right)\mathbf{u}_S + \left(\frac{y_S - y_s}{y_S - y_P}\right)\mathbf{u}_P \right]$ Notice that all the coefficients are scalars. Now, these can be assembled in the form shown below. The rest of the discussion in this thread follows from here. $a_P \mathbf{u}_P + \sum_{nb}{a_{nb} \mathbf{u}_{nb}}$ Hope this helps. Cheers, USV ehsanakrami, rod_petrone, Kummi and 6 others like this.

November 9, 2016, 11:55		#9
Duo New Member DUO ZHANG Join Date: Sep 2016 Posts: 2 Rep Power: 0	thanks, very helpful!

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
What's meaning of UDF FUNCTION	zhaoxinyu	Fluent UDF and Scheme Programming	0	March 31, 2010 09:04
wall y+: physic meaning	sarav	Main CFD Forum	3	November 16, 2009 05:45
want to know meaning	Sangamesh	Siemens	0	May 15, 2007 06:15
What's the meaning of "combustion scalar"and....	cfdbeginner	CFX	0	November 27, 2003 10:02
meaning	id	Main CFD Forum	0	September 18, 2003 05:01