Welcome to TQMesh Discussions!

[FloSewn]

👋 Welcome!

TQMesh is only a small project and still has some bugs and a lack of documentation.
Feel free to ask any questions or share ideas on how to improve it!

[aabramo]

Dear FloSewn,

VERY GOOD PIECE OF WORK!
We (me and colleagues at University of Udine, Italy) are developing a simulator for carrier transport in graphene, and we bumped into your code, which we adopted.

We face a problem.
We are trying to create a device which has an external layer inside which there is another layer and finaly an additional, internal one. However, we wil be interested in adding a FURTHER internal layer.
To summarize: so far we have three nested layers, but in the next future we would like to have a further fourth inside the third one. Unfotunately, we realized that in the 4 nested case the last meshing fails in defining the boundary of the innermost shape, then in meshing it. (at least, this is our understanding)

Are you able to give us a hint on how fix the code? Or - let me dare abusing of your kind support - we can provide you the ".para" file we are using and you may want to check the problem yourself.

I'd be glad if you could reply me with some support.
To this purpose, it would be ok if you replied here or, possibilbly, if you could kindly send me e.mail to: antonio.abramo@uniud.it

Thanks for whatever support you'll be willing to give us.

Cheers,
Antonio Abramo

[HANZhengye]

What's the algrithm of the following code? Would you please offer some reference about this?

inline double evaluate(const Vec2d& xy,
const Domain& domain) const
{
// The underlying size function
const double h_fun = f_(xy);

if ( h_fun <= 0.0 )
  TERMINATE("SizeFunction::evaluate(): Encountered invalid value (<=0).");

double h = h_fun;

// Gather distance contribution of each boundary vertices
for ( const auto& boundary : domain )
  for ( const auto& edge : boundary.get()->edges() )
  {
    const Vec2d&  v1_xy = edge->v1().xy();
    const Vec2d&  v2_xy = edge->v2().xy();

    const double el = edge->length();
    const double r = MAX(h_fun/el, el/h_fun);

    const double d1_sqr = (xy - v1_xy).norm_sqr();
    const double d2_sqr = (xy - v2_xy).norm_sqr();

    const double s1 = (edge->v1().size_range() <= 0.0) 
                    ? el : edge->v1().size_range();
    const double s2 = (edge->v2().size_range() <= 0.0) 
                    ? el : edge->v2().size_range();

    const double s1_inv = 1.0 / (r* s1);
    const double s2_inv = 1.0 / (r* s2);

    const double z1 = exp(-d1_sqr * s1_inv * s1_inv);
    const double z2 = exp(-d2_sqr * s2_inv * s2_inv);

    const double h1 = (edge->v1().mesh_size() <= 0.0)
                    ? h_fun : edge->v1().mesh_size();

    const double h2 = (edge->v2().mesh_size() <= 0.0)
                    ? h_fun : edge->v2().mesh_size();

    const double hv1 = z1*MIN(h1,el) + (1.0-z1)*h_fun;
    const double hv2 = z2*MIN(h2,el) + (1.0-z2)*h_fun;

    h = MIN(h, hv1);
    h = MIN(h, hv2);
  }

// Gather distance contribution of fixed vertices
for ( auto& vertex : domain.fixed_vertices() )
{
  if (vertex->mesh_size() <= 0.0) 
    continue;

  const double d_sqr = (xy - vertex->xy()).norm_sqr();
  const double s_inv = (vertex->size_range() <= 0.0) 
                     ? 1.0/h_fun : 1.0/vertex->size_range();

  const double z = exp(-d_sqr * s_inv * s_inv);
  h = MIN(h, z*vertex->mesh_size() + (1.0-z)*h_fun);
}

return h;

} // evaluate()