Creating rotational symmetric EB around z-axis

[sjlienge]

Dear PeleLMex community, I would like to consult your expertise regarding the setup of user-defined EB. I succeeded to include the EBUserDefined.H file into my setup following the "ChallengeProblem" and want to implement a custom geometry (shown in image) at the inlet which is just a trapezoid rotated around the z-axis.

For the domain I set the following size in the input file:

geometry.prob_lo     = -0.05 -0.05 0.00        # x_lo y_lo (z_lo)
geometry.prob_hi     =  0.05  0.05 0.10        # x_hi y_hi (z_hi)

To generate the shape, in the EBUserDefined.H I want to generate a line from the coordinates x,y,z (0.03,0.0,0.055) and (0.04,0.0,0.05) and then rotate it around the z-axis as in the ChallengeProblem using this setup:

EBUserDefined(
  const Geometry& geom,
  const int required_coarsening_level,
  const int max_coarsening_level)
{
  // ParmParse your geometry parameters
  amrex::ParmParse pp("eb2");

  // Build geometry pieces using EB2::* methods
  EB2::SplineIF Bluff;
  std::vector<amrex::RealVect> lnpts;
  amrex::RealVect p;

  // Define starting and ending point of vector (x,y,z) to describe the line to be rotated
  p = amrex::RealVect(D_DECL(
    0.04, 0.0, 0.05));
  lnpts.push_back(p);
  p = amrex::RealVect(D_DECL(
    0.03, 0.0, 0.055));
  lnpts.push_back(p);
  Bluff.addLineElement(lnpts);
  lnpts.clear();
  // Rotate all objects
  auto revolveBluff = EB2::lathe(Bluff);
  // Build your geometry shop using EB2::makeShop
  auto gshop = EB2::makeShop(revolveBluff);
  // Build geom using EB2::Build
  EB2::Build(gshop, geom, required_coarsening_level, max_coarsening_level);
}

However, this does not work. In fact, the coordinates for p do not appear to be of x,y,z order and they also do not seem to have the same units as in the input file. The whole domain is covered using the coordinates like this. If I change the coordinates to (0.04,0.05,0.0) and (0.03,0.055,0.0) it generates a geomtry, but it is far from what I need. Also how can I make sure only the area below the surface (from lower bound z=0 till surface) is filled?

[sjlienge]

I found now that it is necessary to define the other parts of the cross-section too:

EBUserDefined(
  const Geometry& geom,
  const int required_coarsening_level,
  const int max_coarsening_level)
{
  // ParmParse your geometry parameters
  amrex::ParmParse pp("eb2");

  // Build geometry pieces using EB2::* methods
  EB2::SplineIF Bluff;
  std::vector<amrex::RealVect> lnpts;
  amrex::RealVect p;

  // Define starting and ending point of vector (x,y,z) to describe the line to be rotated

  p = amrex::RealVect(D_DECL(
    0.04, 0.0, 0.00));
  lnpts.push_back(p);
  p = amrex::RealVect(D_DECL(
    0.04, 0.0, 0.05));
  lnpts.push_back(p);
  Bluff.addLineElement(lnpts);
  lnpts.clear();

  p = amrex::RealVect(D_DECL(
    0.04, 0.0, 0.05));
  lnpts.push_back(p);
  p = amrex::RealVect(D_DECL(
    0.03, 0.0, 0.055));
  lnpts.push_back(p);
  Bluff.addLineElement(lnpts);
  lnpts.clear();

  p = amrex::RealVect(D_DECL(
    0.03, 0.0, 0.055));
  lnpts.push_back(p);
  p = amrex::RealVect(D_DECL(
    0.03, 0.0, 0.0));
  lnpts.push_back(p);
  Bluff.addLineElement(lnpts);
  lnpts.clear();

  // Rotate all objects
  auto revolveBluff = EB2::lathe(Bluff);
  // Build your geometry shop using EB2::makeShop
  auto gshop = EB2::makeShop(revolveBluff);
  // Build geom using EB2::Build
  EB2::Build(gshop, geom, required_coarsening_level, max_coarsening_level);
}

This works for a rectangular cross-section but now I receive an error:

amrex::Abort::4::amrex::EB2::build_faces: more than 2 cuts not supported !!!
SIGABRT

[drummerdoc]

I'll have to stare at this a bit, but the error is something that we unfortunately expect to plague any "trailing edge" type geometry (think sharp backside of an airfoil). If the body is thin enough and unfortunately placed, and the background uniform grid is coarse enough, there will be a place with the body cuts through a cell in a way that generates two cut cells at the same index. That is not currently allows in AMReX. The typical solution to that problem is to move the geometry slightly to avoid that problem. However, yours will definitely be tricky, since it's is circular. It will take some experimentation and clever graphics to properly diagnose the issue and plan around it. We'll try to experiment here to find the problem and solution.

That said, I'm a little confused by your picture, because it looks as if the geometry is really a triangle in, say, the XZ plane that is rotated around the z-axis, since it appears that one surface (the "inside surface") is parallel to the z axis, one is in the XY plane, and one is angled. Like I said, I need to stare at this a bit to really understand, but wanted to at least touch on the multicell problem.

[sjlienge]

Thank you so much for the quick response :) I apologize for the confusing picture. I was also adjusting the geometry while I was testing the EB generation and didnt update the graphics, my bad! I have updated the picture:

The two-cell issue makes sense and I admit the current resolution is also very coarse. I will do some experiments on adjusting the mesh and moving the geometry to see if it avoids the error message.

[drummerdoc]

I figured that part out but got sidetracked when updating my comment! I have an idea for visualizing the intersections with the grid that should probably help debugging these things. I'll reach out if I want you help me test it.