DtNmap2Dvisc.hh

Go to the documentation of this file.
1 
9 #ifndef dtn_map2d_visc_hh
10 #define dtn_map2d_visc_hh
11 
12 #include "basics/exceptions.hh"
13 #include "toolbox/sequence.hh"
14 #include "space/integral.hh"
15 #include "formula/elementFormulaContainer.hh"
16 #include "formula/formulas2D.hh"
17 #include "formula/bessel.hh"
18 #include "formula/frmE_product.hh"
19 #include "function/vector.hh"
20 #include "operator/sparseMatrix.hh"
21 #include "hp1D/linearForm.hh"
22 
23 namespace concepts {
24 
25 
26  // ***************************************************** lambda_j_fast *
27 
30  Cmplx lambda_j_fast(const Real L,
31  const Real omega,
32  const Real c,
33  const Real rho0,
34  const Real nu,
35  const uint j)
36  {
37  Cmplx lambda;
38  Real a,b;
39  a = (2.0*j*M_PI/L)*(2.0*j*M_PI/L);
40  b = -omega*rho0/nu;
41 #if __cplusplus >= 201103L
42  lambda.real(-1.0/sqrt(2.0)*sqrt(sqrt(a*a+b*b)+a));
43  lambda.imag(1.0/sqrt(2.0)*sqrt(sqrt(a*a+b*b)-a));
44 #else
45  lambda.real() = -1.0/sqrt(2.0)*sqrt(sqrt(a*a+b*b)+a);
46  lambda.imag() = 1.0/sqrt(2.0)*sqrt(sqrt(a*a+b*b)-a);
47 #endif
48  return lambda;
49  }
50 
51  // ***************************************************** lambda_j_slow *
52 
55  Cmplx lambda_j_slow(const Real L,
56  const Real omega,
57  const Real c,
58  const Real rho0,
59  const Real nu,
60  const uint j)
61  {
62  Cmplx lambda;
63  Real a,b;
64  a = (2.0*j*M_PI/L)*(2.0*j*M_PI/L)
65  - (c*c*omega*omega)/(c*c*c*c+omega*omega*nu*nu/(rho0*rho0));
66  b = -(nu/rho0*omega*omega*omega)/(c*c*c*c+omega*omega*nu*nu/(rho0*rho0));
67 #if __cplusplus >= 201103L
68  lambda.real(-1.0/sqrt(2.0)*sqrt(sqrt(a*a+b*b)+a));
69  lambda.imag(1.0/sqrt(2.0)*sqrt(sqrt(a*a+b*b)-a));
70 #else
71  lambda.real() = -1.0/sqrt(2.0)*sqrt(sqrt(a*a+b*b)+a);
72  lambda.imag() = 1.0/sqrt(2.0)*sqrt(sqrt(a*a+b*b)-a);
73 #endif
74  return lambda;
75  }
76 
77  // ***************************************************** g_fast *
80  Cmplx g_fast(const Real L,
81  const Real omega,
82  const Real c,
83  const Real rho0,
84  const Real nu,
85  const uint j)
86  {
87  Cmplx g;
88 #if __cplusplus >= 201103L
89  g.real(-c*c*rho0*(2.0*j*M_PI/L)*(2.0*j*M_PI/L));
90  g.imag(0.);
91 #else
92  g.real() = -c*c*rho0*(2.0*j*M_PI/L)*(2.0*j*M_PI/L);
93  g.imag() = 0.;
94 #endif
95  return g;
96  }
97 
98  // ***************************************************** g_slow *
101  Cmplx g_slow(const Real L,
102  const Real omega,
103  const Real c,
104  const Real rho0,
105  const Real nu,
106  const uint j)
107  {
108  Cmplx g;
109 #if __cplusplus >= 201103L
110  g.real(omega*omega*rho0-c*c*rho0*(2.0*j*M_PI/L)*(2.0*j*M_PI/L)
111  - (omega*omega*omega*omega*nu*nu/rho0)/(c*c*c*c+nu*nu*omega*omega/(rho0*rho0)));
112  g.imag(omega*omega*omega*nu*c*c/(c*c*c*c+nu*nu*omega*omega/(rho0*rho0)));
113 #else
114  g.real() = omega*omega*rho0-c*c*rho0*(2.0*j*M_PI/L)*(2.0*j*M_PI/L)
115  - (omega*omega*omega*omega*nu*nu/rho0)/(c*c*c*c+nu*nu*omega*omega/(rho0*rho0));
116  g.imag() = omega*omega*omega*nu*c*c/(c*c*c*c+nu*nu*omega*omega/(rho0*rho0));
117 #endif
118  return g;
119  }
120 
121 
122 
123  Sequence<Cmplx> getNSDtNCoeff_Straight2D_partDn(const Real L,
124  const Real omega,
125  const Real c,
126  const Real rho0,
127  const Real nu,
128  uint N=0)
129  {
130 
131  Sequence<int> nSeq = makeRangeSequence(N+1);
132  Sequence<Cmplx> DtNCoeff(N+1);
133  DtNCoeff[0] = lambda_j_slow(L,omega,c,rho0,nu,0);
134  for(uint j=1;j <= N; j++)
135  {
136  Cmplx ljf = lambda_j_fast(L,omega,c,rho0,nu,j);
137  Cmplx ljs = lambda_j_slow(L,omega,c,rho0,nu,j);
138  Cmplx gjf = g_fast(L,omega,c,rho0,nu,j);
139  Cmplx gjs = g_slow(L,omega,c,rho0,nu,j);
140  DtNCoeff[j] = ljf*ljs*(gjf-gjs)
141  /(gjf*ljs-gjs*ljf);
142  }
143  return DtNCoeff;
144  }
145 
146  Sequence<Cmplx> getNSDtNCoeff_Straight2D_partDt(const Real L,
147  const Real omega,
148  const Real c,
149  const Real rho0,
150  const Real nu,
151  uint N=1)
152  {
153 
154  Sequence<int> nSeq = makeRangeSequence(N);
155  Sequence<Cmplx> DtNCoeff(N);
156  for(uint j=0;j < N; j++)
157  {
158  Cmplx ljf = lambda_j_fast(L,omega,c,rho0,nu,j+1);
159  Cmplx ljs = sqrt(nu)*lambda_j_slow(L,omega,c,rho0,nu,j+1);
160  Cmplx gjf = g_fast(L,omega,c,rho0,nu,j+1);
161  Cmplx gjs = g_slow(L,omega,c,rho0,nu,j+1);
162  DtNCoeff[j] = 2.0 * (j+1) * M_PI / L
163  - 1.0 / (c*c*rho0 * 2.0 * (j+1) * M_PI / L)
164  * ((ljf-ljs)*gjf*gjs) / (gjf*ljs-gjs*ljf);
165  }
166  return DtNCoeff;
167  }
168 
169  Sequence<Cmplx> getNSDtNCoeff_Straight2D_partDttilde(const Real L,
170  const Real omega,
171  const Real c,
172  const Real rho0,
173  const Real nu,
174  uint N=1)
175  {
176 
177  Sequence<int> nSeq = makeRangeSequence(N);
178  Sequence<Cmplx> DtNCoeff(N);
179  for(uint j=0;j < N; j++)
180  {
181  Cmplx ljf = lambda_j_fast(L,omega,c,rho0,nu,j+1);
182  Cmplx ljs = sqrt(nu)*lambda_j_slow(L,omega,c,rho0,nu,j+1);
183  Cmplx gjf = g_fast(L,omega,c,rho0,nu,j+1);
184  Cmplx gjs = g_slow(L,omega,c,rho0,nu,j+1);
185  DtNCoeff[j] = - 1.0 / (c*c*rho0 *2.0 * (j+1) * M_PI / L)
186  * ((ljf-ljs)*gjf*gjs) / (gjf*ljs-gjs*ljf);
187  }
188  return DtNCoeff;
189  }
190 
191  Sequence<Cmplx> getNSDtNCoeff_Straight2D_partRn(const Real L,
192  const Real omega,
193  const Real c,
194  const Real rho0,
195  const Real nu,
196  uint N=1)
197  {
198 
199  Sequence<int> nSeq = makeRangeSequence(N);
200  Sequence<Cmplx> DtNCoeff(N);
201  for(uint j=0;j < N; j++)
202  {
203  Cmplx ljf = lambda_j_fast(L,omega,c,rho0,nu,j+1);
204  Cmplx ljs = lambda_j_slow(L,omega,c,rho0,nu,j+1);
205  Cmplx gjf = g_fast(L,omega,c,rho0,nu,j+1);
206  Cmplx gjs = g_slow(L,omega,c,rho0,nu,j+1);
207  DtNCoeff[j] = nu * 2.0 * (j+1) * M_PI / L
208  + c*c*rho0 * 2.0 * (j+1) * M_PI / L
209  * nu *((ljf-ljs)*ljs*ljf) / (gjf*ljs-gjs*ljf);
210  }
211  return DtNCoeff;
212  }
213 
214  Sequence<Cmplx> getNSDtNCoeff_Straight2D_partRntilde(const Real L,
215  const Real omega,
216  const Real c,
217  const Real rho0,
218  const Real nu,
219  uint N=1)
220  {
221 
222  Sequence<int> nSeq = makeRangeSequence(N);
223  Sequence<Cmplx> DtNCoeff(N);
224  for(uint j=0;j < N; j++)
225  {
226  Cmplx ljf = lambda_j_fast(L,omega,c,rho0,nu,j+1);
227  Cmplx ljs = lambda_j_slow(L,omega,c,rho0,nu,j+1);
228  Cmplx gjf = g_fast(L,omega,c,rho0,nu,j+1);
229  Cmplx gjs = g_slow(L,omega,c,rho0,nu,j+1);
230  DtNCoeff[j] = c*c*rho0 * 2.0 * (j+1) * M_PI / L
231  * nu *((ljf-ljs)*ljs*ljf) / (gjf*ljs-gjs*ljf);
232  }
233  return DtNCoeff;
234  }
235 
236  Sequence<Cmplx> getNSDtNCoeff_Straight2D_partRt(const Real L,
237  const Real omega,
238  const Real c,
239  const Real rho0,
240  const Real nu,
241  uint N=0)
242  {
243 
244  Sequence<int> nSeq = makeRangeSequence(N+1);
245  Sequence<Cmplx> DtNCoeff(N+1);
246  DtNCoeff[0] = nu*lambda_j_fast(L,omega,c,rho0,nu,0);
247  for(uint j=1;j <= N; j++)
248  {
249  Cmplx ljf = lambda_j_fast(L,omega,c,rho0,nu,j);
250  Cmplx ljs = lambda_j_slow(L,omega,c,rho0,nu,j);
251  Cmplx gjf = g_fast(L,omega,c,rho0,nu,j);
252  Cmplx gjs = g_slow(L,omega,c,rho0,nu,j);
253  DtNCoeff[j] = nu*(gjf*ljs*ljs-gjs*ljf*ljf)
254  /(gjf*ljs-gjs*ljf);
255  }
256  return DtNCoeff;
257  }
258 
259 } // namespace concepts
260 
261 #endif // dtn_map2d_visc_hh
concepts::getNSDtNCoeff_Straight2D_partDn
Sequence< Cmplx > getNSDtNCoeff_Straight2D_partDn(const Real L, const Real omega, const Real c, const Real rho0, const Real nu, uint N=0)
Definition: DtNmap2Dvisc.hh:123
concepts::getNSDtNCoeff_Straight2D_partRn
Sequence< Cmplx > getNSDtNCoeff_Straight2D_partRn(const Real L, const Real omega, const Real c, const Real rho0, const Real nu, uint N=1)
Definition: DtNmap2Dvisc.hh:191
concepts::getNSDtNCoeff_Straight2D_partRntilde
Sequence< Cmplx > getNSDtNCoeff_Straight2D_partRntilde(const Real L, const Real omega, const Real c, const Real rho0, const Real nu, uint N=1)
Definition: DtNmap2Dvisc.hh:214
M_PI
#define M_PI
Definition: typedefs.hh:13
concepts::getNSDtNCoeff_Straight2D_partRt
Sequence< Cmplx > getNSDtNCoeff_Straight2D_partRt(const Real L, const Real omega, const Real c, const Real rho0, const Real nu, uint N=0)
Definition: DtNmap2Dvisc.hh:236
concepts::Sequence
Sequence with operations, output operator, and method of the particular element types.
Definition: sequence.hh:39
concepts::Cmplx
std::complex< Real > Cmplx
Type for a complex number. It also depends on the setting of Real.
Definition: typedefs.hh:39
concepts::getNSDtNCoeff_Straight2D_partDt
Sequence< Cmplx > getNSDtNCoeff_Straight2D_partDt(const Real L, const Real omega, const Real c, const Real rho0, const Real nu, uint N=1)
Definition: DtNmap2Dvisc.hh:146
concepts::getNSDtNCoeff_Straight2D_partDttilde
Sequence< Cmplx > getNSDtNCoeff_Straight2D_partDttilde(const Real L, const Real omega, const Real c, const Real rho0, const Real nu, uint N=1)
Definition: DtNmap2Dvisc.hh:169
concepts::lambda_j_fast
Cmplx lambda_j_fast(const Real L, const Real omega, const Real c, const Real rho0, const Real nu, const uint j)
Compute the eigenvalues.
Definition: DtNmap2Dvisc.hh:30
concepts::lambda_j_slow
Cmplx lambda_j_slow(const Real L, const Real omega, const Real c, const Real rho0, const Real nu, const uint j)
Compute the eigenvalues.
Definition: DtNmap2Dvisc.hh:55
concepts::g_slow
Cmplx g_slow(const Real L, const Real omega, const Real c, const Real rho0, const Real nu, const uint j)
Compute the coefficients.
Definition: DtNmap2Dvisc.hh:101
concepts::makeRangeSequence
Sequence< int > makeRangeSequence(uint n)
Returns the sequence 0,1,...,n-1.
concepts::g_fast
Cmplx g_fast(const Real L, const Real omega, const Real c, const Real rho0, const Real nu, const uint j)
Compute the coefficients.
Definition: DtNmap2Dvisc.hh:80
frmE_product.hh
concepts::Real
double Real
Type normally used for a floating point number.
Definition: typedefs.hh:36
concepts
Basic namespace for Concepts-2.
Definition: pml_formula.h:16
Page URL: http://wiki.math.ethz.ch/bin/view/Concepts/WebHome
21 August 2020
© 2020 Eidgenössische Technische Hochschule Zürich