LifeV
OneDFSIFluxNonLinear.cpp
Go to the documentation of this file.
1 //@HEADER
2 /*
3 *******************************************************************************
4 
5  Copyright (C) 2004, 2005, 2007 EPFL, Politecnico di Milano, INRIA
6  Copyright (C) 2010 EPFL, Politecnico di Milano, Emory University
7 
8  This file is part of LifeV.
9 
10  LifeV is free software; you can redistribute it and/or modify
11  it under the terms of the GNU Lesser General Public License as published by
12  the Free Software Foundation, either version 3 of the License, or
13  (at your option) any later version.
14 
15  LifeV is distributed in the hope that it will be useful,
16  but WITHOUT ANY WARRANTY; without even the implied warranty of
17  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
18  Lesser General Public License for more details.
19 
20  You should have received a copy of the GNU Lesser General Public License
21  along with LifeV. If not, see <http://www.gnu.org/licenses/>.
22 
23 *******************************************************************************
24 */
25 //@HEADER
26 
27 /*!
28  * @file
29  * @brief File containing a base class for non linear 1D model flux function.
30  *
31  * @version 1.0
32  * @author Vincent Martin
33  *
34  * @version 2.0
35  * @date 15-04-2010
36  * @author Cristiano Malossi <cristiano.malossi@epfl.ch>
37  *
38  * @contributors Simone Rossi <simone.rossi@epfl.ch>, Ricardo Ruiz-Baier <ricardo.ruiz@epfl.ch>
39  * @maintainer Cristiano Malossi <cristiano.malossi@epfl.ch>
40  */
41 
42 #include <lifev/one_d_fsi/solver/OneDFSIFluxNonLinear.hpp>
43 
44 namespace LifeV
45 {
46 // ===================================================
47 // Methods
48 // ===================================================
49 Real
50 OneDFSIFluxNonLinear::flux ( const Real& A, const Real& Q, const ID& row, const UInt& iNode ) const
51 {
52  if ( row == 0 ) // F1
53  {
54  return Q;
55  }
56 
57  if ( row == 1 ) // F2
58  {
59  return ( M_physicsPtr->data()->alpha ( iNode ) * Q * Q / A +
60  M_physicsPtr->data()->beta0 ( iNode ) * M_physicsPtr->data()->beta1 ( iNode ) *
61  M_physicsPtr->data()->area0 ( iNode ) / ( ( M_physicsPtr->data()->beta1 ( iNode ) + 1 ) *
62  M_physicsPtr->data()->densityRho() ) * ( OneDFSI::pow15 ( A / M_physicsPtr->data()->area0 ( iNode ),
63  M_physicsPtr->data()->beta1 ( iNode ) + 1 ) - 1 ) ) *
64  M_physicsPtr->data()->robertsonCorrection();
65  }
66 
67  ERROR_MSG ("The flux function has only 2 components.");
68 
69  return -1.;
70 }
71 
72 Real
73 OneDFSIFluxNonLinear::dFdU ( const Real& A, const Real& Q, const ID& row, const ID& column, const UInt& iNode ) const
74 {
75  if ( row == 0 && column == 0 ) // dF1/dA
76  {
77  return 0.;
78  }
79 
80  if ( row == 0 && column == 1 ) // dF1/dQ
81  {
82  return 1.;
83  }
84 
85  if ( row == 1 && column == 0 ) // dF2/dA
86  {
87  return ( M_physicsPtr->data()->beta0 ( iNode ) *
88  M_physicsPtr->data()->beta1 ( iNode ) /
89  M_physicsPtr->data()->densityRho() * OneDFSI::pow05 ( A / M_physicsPtr->data()->area0 ( iNode ),
90  M_physicsPtr->data()->beta1 ( iNode ) ) -
91  M_physicsPtr->data()->alpha ( iNode ) * Q * Q / A / A ) *
92  M_physicsPtr->data()->robertsonCorrection();
93  }
94  if ( row == 1 && column == 1 ) // dF2/dQ
95  {
96  return M_physicsPtr->data()->robertsonCorrection() * 2 * M_physicsPtr->data()->alpha ( iNode ) * Q / A;
97  }
98 
99  ERROR_MSG ("Flux's differential function has only 4 components.");
100 
101  return -1.;
102 }
103 
104 void
105 OneDFSIFluxNonLinear::eigenValuesEigenVectors ( const Real& A,
106  const Real& Q,
107  container2D_Type& eigenvalues,
108  container2D_Type& leftEigenvector1,
109  container2D_Type& leftEigenvector2,
110  const UInt& iNode ) const
111 {
112 #ifdef HAVE_LIFEV_DEBUG
113  debugStream (6312) << "[OneDFSIModel_Flux_NonLinear]::jacobian_EigenValues_Vectors\n";
114 #endif
115 
116  Real celerity;
117  celerity = std::sqrt ( M_physicsPtr->data()->alpha ( iNode ) * (
118  M_physicsPtr->data()->alpha ( iNode ) - 1) * Q * Q / ( A * A ) +
119  M_physicsPtr->data()->beta0 ( iNode ) *
120  M_physicsPtr->data()->beta1 ( iNode ) /
121  M_physicsPtr->data()->densityRho() * OneDFSI::pow05 ( A / M_physicsPtr->data()->area0 ( iNode ),
122  M_physicsPtr->data()->beta1 ( iNode ) ) );
123 
124  eigenvalues[0] = M_physicsPtr->data()->alpha ( iNode ) * Q / A + celerity;
125  eigenvalues[1] = M_physicsPtr->data()->alpha ( iNode ) * Q / A - celerity;
126 
127  leftEigenvector1[0] = - eigenvalues[1] / A;
128  leftEigenvector1[1] = 1. / A;
129  leftEigenvector2[0] = - eigenvalues[0] / A;
130  leftEigenvector2[1] = 1. / A;
131 }
132 
133 void
134 OneDFSIFluxNonLinear::deltaEigenValuesEigenVectors ( const Real& A,
135  const Real& Q,
136  container2D_Type& deltaEigenvalues,
137  container2D_Type& deltaLeftEigenvector1,
138  container2D_Type& deltaLeftEigenvector2,
139  const UInt& iNode ) const
140 {
141  Real deltaCelerity;
142 
143  Real AoverA0 ( A / M_physicsPtr->data()->area0 ( iNode ) );
144  Real C ( OneDFSI::pow05 ( AoverA0, M_physicsPtr->data()->beta1 ( iNode ) ) / M_physicsPtr->data()->densityRho() );
145 
146  deltaCelerity = 0.5 / std::sqrt ( M_physicsPtr->data()->alpha ( iNode ) * (
147  M_physicsPtr->data()->alpha ( iNode ) - 1) * Q * Q / ( A * A ) +
148  M_physicsPtr->data()->beta0 ( iNode ) *
149  M_physicsPtr->data()->beta1 ( iNode ) * C ) * ( C * (
150  M_physicsPtr->data()->beta1 ( iNode ) *
151  M_physicsPtr->data()->dBeta0dz ( iNode ) -
152  M_physicsPtr->data()->beta0 ( iNode ) *
153  M_physicsPtr->data()->beta1 ( iNode ) *
154  M_physicsPtr->data()->beta1 ( iNode ) /
155  M_physicsPtr->data()->area0 ( iNode ) *
156  M_physicsPtr->data()->dArea0dz ( iNode ) +
157  M_physicsPtr->data()->beta0 ( iNode ) * ( 1 +
158  M_physicsPtr->data()->beta0 ( iNode ) * std::log ( AoverA0 ) ) *
159  M_physicsPtr->data()->dBeta1dz ( iNode ) ) + ( 2 *
160  M_physicsPtr->data()->alpha ( iNode ) - 1 ) * Q * Q / ( A * A ) *
161  M_physicsPtr->data()->dAlphadz ( iNode ) );
162 
163  deltaEigenvalues[0] = M_physicsPtr->data()->dAlphadz ( iNode ) * Q / A + deltaCelerity;
164  deltaEigenvalues[1] = M_physicsPtr->data()->dAlphadz ( iNode ) * Q / A - deltaCelerity;
165 
166  deltaLeftEigenvector1[0] = - deltaEigenvalues[1] / A;
167  deltaLeftEigenvector1[1] = 0;
168  deltaLeftEigenvector2[0] = - deltaEigenvalues[0] / A;
169  deltaLeftEigenvector2[1] = 0;
170 }
171 
172 }
LifeV::OneDFSIFluxNonLinear::eigenValuesEigenVectors
void eigenValuesEigenVectors(const Real &A, const Real &Q, container2D_Type &eigenvalues, container2D_Type &leftEigenvector1, container2D_Type &leftEigenvector2, const UInt &iNode) const
Eigenvalues and eigenvectors of the Jacobian matrix.
Definition: OneDFSIFluxNonLinear.cpp:105
LifeV::OneDFSIFluxNonLinear::dFdU
Real dFdU(const Real &A, const Real &Q, const ID &row, const ID &column, const UInt &iNode) const
Evaluate the derivative of the flux term.
Definition: OneDFSIFluxNonLinear.cpp:73
ETCurrentFE::updateInverseJacobian
void updateInverseJacobian(const UInt &iQuadPt)
Definition: ETCurrentFE.cpp:405
LifeV::OneDFSIFluxNonLinear::deltaEigenValuesEigenVectors
void deltaEigenValuesEigenVectors(const Real &A, const Real &Q, container2D_Type &deltaEigenvalues, container2D_Type &deltaLeftEigenvector1, container2D_Type &deltaLeftEigenvector2, const UInt &iNode) const
Derivatives of the eigenvalues and eigenvectors of the derivative of the Jacobian matrix...
Definition: OneDFSIFluxNonLinear.cpp:134
ERROR_MSG
#define ERROR_MSG(A)
Definition: LifeAssert.hpp:69
LifeV::ID
uint32_type ID
IDs.
Definition: LifeV.hpp:194
LifeV::OneDFSIFluxNonLinear
OneDFSIFluxNonLinear - Class containing the non-linear flux term of the 1D hyperbolic problem...
Definition: OneDFSIFluxNonLinear.hpp:108
LifeV::Real
double Real
Generic real data.
Definition: LifeV.hpp:175
LifeV::OneDFSIFluxNonLinear::flux
Real flux(const Real &A, const Real &Q, const ID &row, const UInt &iNode) const
Evaluate the flux term.
Definition: OneDFSIFluxNonLinear.cpp:50
LifeV::UInt
uint32_type UInt
generic unsigned integer (used mainly for addressing)
Definition: LifeV.hpp:191