Fuerza de arrastre y sustentación que actúa sobre un perfil aerodinámico

Esta duda se puede abordar de variadas maneras, pero nosotros te damos la respuesta más completa en nuestra opinión.

Solución:

Cuando ejecuto tu código obtengo un FindRoot mensaje de advertencia:

Lo que me hace sospechar de la calidad del resultado. Si asumimos que el resultado es correcto, podemos acelerar la integración usando el FEM también para eso. Creamos una malla de elemento límite de la lámina:

bmeshFoil = 
  ToBoundaryMesh["Coordinates" -> coords[[5 ;; nn]], 
   "BoundaryElements" -> LineElement[
      Partition[Range[Length[coords[[5 ;; nn]]]], 2, 1, 1]]];

E integre a lo largo del límite:

fdrag, flift = 
 NIntegrate[force[x, y], x, y [Element] bmeshFoil, 
   AccuracyGoal -> 3, PrecisionGoal -> 3] // AbsoluteTiming

(* 0.702661, 0.209457, 1.34502 *)

Aquí hay un no parcialNIntegrate respuesta que aún necesita trabajo pero que podría darle algunas ideas sobre cómo proceder.

Extendí el dominio para que me fuera más fácil elegir segmentos de línea relacionados con el perfil aerodinámico.

x1 = -2; x2 = 3; y1 = -1.5; y2 = 1.5;(*domain dimensions*)

Luego seguí este ejemplo de la documentación para tomar las normales en el punto medio del segmento de línea y la longitud de cada segmento:

bn = bmesh["BoundaryNormals"];
mean = Mean /@ GetElementCoordinates[bmesh["Coordinates"], #] & /@ 
   ElementIncidents[bmesh["BoundaryElements"]];
dist = EuclideanDistance @@@ 
     GetElementCoordinates[bmesh["Coordinates"], #] & /@ 
   ElementIncidents[bmesh["BoundaryElements"]];
ids = [email protected]






   Position[
    Flatten[mean, 1], _?(EuclideanDistance[#, 0, 0] < 1.1 &), 1];
foilbn = bn[[1, ids]];
foilbnplt = ArrayReshape[foilbn, 1~Join~(foilbn // Dimensions)];
foildist = dist[[1, ids]];
foildistplt = 
  ArrayReshape[foildist, 1~Join~(foildist // Dimensions)];
foilmean = mean[[1, ids]];
foilmeanplt = 
  ArrayReshape[foilmean, 1~Join~(foilmean // Dimensions)];
Show[bmesh["Wireframe"], 
 Graphics[MapThread[
   Arrow[#1, #2] &, Join @@ foilmeanplt, 
    Join @@ (foilbnplt/5 + foilmeanplt)]]]

Parece que capturamos todas las normales asociadas con el perfil aerodinámico. Tienes muchas normales, así que creo que una suma ponderada debería ser una aproximación decente a la integral.

Luego, creé una función que toma una suma ponderada de fuerzas. Es rápido pero necesita algo de trabajo y validación, pero este método es similar a lo que se hace con otros códigos.

ClearAll[fluidStress]
fluidStress[uif_InterpolatingFunction, vif_InterpolatingFunction, 
   pif_InterpolatingFunction, mu_, rho_, bn_, dist_, mean_] := 
 Block[dd, df, mesh, coords, dv, press, fx, fy, wfx, wfy, nx, ny, ux,
    uy, vx, vy, 
  dd = Outer[(D[#1[x, y], #2]) &, uif, vif, x, y];
  df = Table[Function[x, y, Evaluate[dd[[i, j]]]], i, 2, j, 2];
  (*the coordinates from the foil*)
  coords = mean;
  dv = Table[df[[i, j]] @@@ coords, i, 2, j, 2];
  ux = dv[[1, 1]];
  uy = dv[[1, 2]];
  vx = dv[[2, 1]];
  vy = dv[[2, 2]];
  nx = bn[[All, 1]];
  ny = bn[[All, 2]];
  press = pif[#1, #2] & @@@ coords;
  fx = -nx*press + mu*(-2*nx*ux - ny*(uy + vx));
  fy = -ny*press + mu*(-nx*(vx + uy) - 2*ny*vy);
  wfx = dist*fx ;
  wfy = dist*fy; 
  Total /@ wfx, wfy
  ]
AbsoluteTiming[fdrag, flift = 
  fluidStress[xVel, yVel, pressure, 10^-3, 1, foilbn, foildist, 
   foilmean]]
(* 0.364506, 0.00244262, 0.158859 *)

Sección de Reseñas y Valoraciones

Eres capaz de añadir valor a nuestra información añadiendo tu experiencia en las acotaciones.

---

---

Utiliza Nuestro Buscador