update example

examples/helmholtzflow/Nacelle/main.cpp

@@ -28,7 +28,7 @@ int main(int argc, char **argv)
 {
   GmshDdm gmshDdm(argc, argv);
 
-  int ElemOrder = 4; // mesh order
+  int ElemOrder = 2; // mesh order
   double h_near = 0.03;
   double h_far = 0.04;
   double Lx = 3.;
@@ -125,7 +125,7 @@ int main(int argc, char **argv)
   gmsh::model::mesh::setOrder(ElemOrder);
   gmsh::write("nacelle_post.msh");
 
-  bool saveFlow = 0;
+  bool saveFlow = 1;
   // Compute mean flow properties, Poisson problem
   // Mean flow parameters
   double gamma0 = 1.4;
@@ -190,10 +190,7 @@ int main(int argc, char **argv)
   rho0 = rho_inf * pow( (1 + (gamma0-1)/2.* (v_inf*v_inf - v*v) / (c_inf*c_inf) ), 1/(gamma0-1) );
 
   if (saveFlow) {
-    save(rho0, omegaFlow, "rho0");
-    save(c0, omegaFlow, "c0");
-    save(M0x, omegaFlow, "Mx");
-    save(M0y, omegaFlow, "My");
+    save(sqrt(M0x*M0x+M0y*M0y), omegaFlow, "Mach"+std::to_string(ElemOrder));
   }
 
   // DDM problem
@@ -223,7 +220,7 @@ int main(int argc, char **argv)
   getInputMode(m, n, w/cf, Mf, kx, ky, psi, dy_psi, dz_psi, 2);
   msg::info << "propagative wavenumber kx = " << kx << ", transverse wavenumber ky = " << ky << msg::endl;
 
-  float pointsByWl = ( 2*pi*FEMorder*c_inf*( 1 - abs(Mf) ) ) / (h_far*w);
+  float pointsByWl = ( 2*pi*FEMorder*cf*( 1 - abs(Mf) ) ) / (h_near*w);
   msg::info << " - FEM basis order = " << FEMorder << "" << msg::endl;
   msg::info << " - Approx. dofs by wavelength at fan Face = " << pointsByWl << msg::endl;
 
@@ -311,7 +308,7 @@ int main(int argc, char **argv)
     // convected Helmholz weak form
     formulation(i).integral(rho0S * vector< std::complex< double > >(1-Mx*Mx, -Mx*My, 0.) * grad(dof(u(i))), vector< std::complex< double > >(1., 0., 0.) * grad(tf(u(i))), omega_phy(i), gauss);
     formulation(i).integral(rho0S * vector< std::complex< double > >(-Mx*My, 1-My*My, 0.) * grad(dof(u(i))), vector< std::complex< double > >(0., 1., 0.) * grad(tf(u(i))), omega_phy(i), gauss);
-    formulation(i).integral(rho0S * dof(u(i)), vector< std::complex< double > >(-im*k0*Mx,-im*k0*My, 0.) * grad(tf(u(i))), omega_phy(i), gauss, term::ProductType::Scalar);
+    formulation(i).integral(-im * k0 * rho0S * dof(u(i)), vector< std::complex< double > >(Mx, My, 0.) * grad(tf(u(i))), omega_phy(i), gauss);
     formulation(i).integral(rho0S * vector< std::complex< double > >(im*k0*Mx, im*k0*My, 0.) * grad(dof(u(i))), tf(u(i)), omega_phy(i), gauss);
     formulation(i).integral(-rho0S * k0 * k0 * dof(u(i)), tf(u(i)), omega_phy(i), gauss);
     formulation(i).integral(+ rho0S *(m*m) / ( y< std::complex< double > >()*y< std::complex< double > >() ) * dof(u(i)), tf(u(i)), omega_phy(i), gauss);