maia::geom Namespace Reference

Functions
template<typename T , MInt nDim>
MBool	doBoxesOverlap (const T *const a, const T *const b)
template<typename T , MInt nDim>
MBool	isBoxInsideBox (const T *const a, const T *const b)
template<typename T , MInt nDim>
MBool	isPointInsideBox (const T *const p, const T *const b)
template<MInt nDim>
MBool	isPointInsideSphere (const MFloat *const p, const MFloat *const c, const MFloat R)
template<MInt nDim>
MBool	isPointInsideRingSegment (const MFloat *const p, const MFloat *const c, const MFloat *const R, const MFloat *const phi, const MInt axis)
template<MInt nDim>
MBool	isBoxInsideSphere (const MFloat *const p, const MFloat *const c, const MFloat R)
template<typename T , MInt nDim>
MBool	isSphereInsideBox (const T *const c, const T R, const T *const p)
template<MInt nDim>
MBool	doBoxAndSphereOverlap (const MFloat *const b, const MFloat *const c, const MFloat R)
template<MInt nDim>
MBool	doesLinePenetrateBox (const MFloat *const line, const MFloat *const bbox)
template<MInt nDim>
MBool	doesLinePenetrateSphere (const MFloat *const line, const MFloat *const center, const MFloat radius)

Function Documentation

doBoxAndSphereOverlap()

template<MInt nDim>

MBool maia::geom::doBoxAndSphereOverlap	(	const MFloat *const	b,
		const MFloat *const	c,
		const MFloat	R
	)

Definition at line 123 of file geometryanalytic.h.

                                                                                          {
  MFloat cp[nDim];
  const MInt cornerCode[8][3] = {{0, 0, 0},    {nDim, 0, 0},    {0, nDim, 0},    {nDim, nDim, 0},
                                 {0, 0, nDim}, {nDim, 0, nDim}, {0, nDim, nDim}, {nDim, nDim, nDim}};
 
  // easy part
  // is sphere center inside box?
  if(maia::geom::isPointInsideBox<MFloat, nDim>(c, b)) {
    return true;
  }
 
  // is a box corner in sphere?
  for(MInt corner = 0; corner < IPOW2(nDim); corner++) {
    for(MInt dim = 0; dim < nDim; dim++) {
      cp[dim] = b[dim + cornerCode[corner][dim]];
    }
    if(maia::geom::isPointInsideSphere<nDim>(cp, c, R)) {
      return true;
    }
  }
 
  // 'tricky part' reducing the dimensions
  // what region relative to the cube is the sphere center located in
  MInt insideCnt = 0;
  MFloat p[nDim];
  MFloat lc[nDim];
  for(MInt dim = 0; dim < nDim; dim++) {
    if(c[dim] > b[dim] && c[dim] < b[dim + nDim]) {
      p[dim] = F0;
      lc[dim] = F0;
      insideCnt++;
    } else {
      lc[dim] = c[dim];
      p[dim] = (c[dim] <= b[dim]) ? b[dim] : b[dim + nDim];
    }
  }
 
  // the closest point is a corner, we already checked them
  if(insideCnt == 0) {
    return false;
  }
 
  // the cosest point is found checking the outside dimensions
  return maia::geom::isPointInsideSphere<nDim>(p, lc, R);
}

doBoxesOverlap()

template<typename T , MInt nDim>

MBool maia::geom::doBoxesOverlap	(	const T *const	a,
		const T *const	b
	)

Definition at line 19 of file geometryanalytic.h.

                                                         {
  for(MInt dim = 0; dim < nDim; dim++) {
    if((a[dim] > b[dim + nDim]) || b[dim] > a[dim + nDim]) {
      return false;
    }
  }
  return true;
}

doesLinePenetrateBox()

template<MInt nDim>

MBool maia::geom::doesLinePenetrateBox	(	const MFloat *const	line,
		const MFloat *const	bbox
	)

Definition at line 170 of file geometryanalytic.h.

                                                                               {
  IF_CONSTEXPR(nDim > 3 || nDim < 2) { mTerm(1, AT_, "neither 2 nor 3 dimensional"); }
 
  // create a line bounding box
  MFloat bbox_line[2 * nDim];
  for(MInt dim = 0; dim < nDim; dim++) {
    if(line[dim] < line[dim + nDim]) {
      bbox_line[dim] = line[dim];
      bbox_line[dim + nDim] = line[dim + nDim];
    } else {
      bbox_line[dim] = line[dim + nDim];
      bbox_line[dim + nDim] = line[dim];
    }
  }
 
  // check bounding boxes for overlap
  // a separation axis can be found in the cartesian axis
  if(!maia::geom::doBoxesOverlap<MFloat, nDim>(bbox_line, bbox)) {
    return false;
  }
 
  // translate box to center and scale to unit cube
  MFloat line_c[2 * nDim];
  for(MInt dim = 0; dim < nDim; dim++) {
    const MFloat halfLength = (bbox[dim + nDim] - bbox[dim]) / F2;
    line_c[dim] = (line[dim] - bbox[dim]) / halfLength - 0.5;
    line_c[dim + nDim] = (line[dim + nDim] - bbox[dim]) / halfLength - 0.5;
  }
 
  MFloat line_span[3] = {F0, F0, F0};
  MFloat lineLength = F0;
  for(MInt dim = 0; dim < nDim; dim++) {
    line_span[dim] = line_c[dim + nDim] - line_c[dim];
    lineLength += POW2(line_span[dim]);
  }
  lineLength = std::sqrt(lineLength);
  for(MInt dim = 0; dim < nDim; dim++) {
    line_span[dim] /= lineLength;
  }
 
  const MFloat radius = nDim * F1B4;
  for(MInt e = 2; e > ((nDim == 2) ? 1 : -1); e--) {
    // generate separation axis
    MFloat axis[3];
    MFloat e_vec[3] = {F0, F0, F0};
    e_vec[e] = F1;
    maia::math::cross(line_span, e_vec, axis);
 
    // project one of the points on the axis
    MFloat line_projection = F0;
    for(MInt dim = 0; dim < nDim; dim++) {
      line_projection += axis[dim] * line_c[dim];
    }
 
    // larger equal radius means no penetration but maybe a contactpoint
    if(POW2(line_projection) >= radius) {
      return false;
    }
 
    // project the corner of the cube ...
    // which corner to project
    MFloat corner_projection = F0;
    for(MInt dim = 0; dim < nDim; dim++) {
      corner_projection += axis[dim] * ((axis[dim] > 0) ? 0.5 : -0.5);
    }
 
    // larger equal radius means no penetration but maybe a contactpoint
    if(std::fabs(line_projection) >= std::fabs(corner_projection)) {
      return false;
    }
  }
 
  // if no separation axis can be found, there has to be a penetration
  return true;
}

doesLinePenetrateSphere()

template<MInt nDim>

MBool maia::geom::doesLinePenetrateSphere	(	const MFloat *const	line,
		const MFloat *const	center,
		const MFloat	radius
	)

Definition at line 247 of file geometryanalytic.h.

                                                                                                         {
  const MBool p1_inside = maia::geom::isPointInsideSphere<nDim>(line, center, radius);
  const MBool p2_inside = maia::geom::isPointInsideSphere<nDim>(line + nDim, center, radius);
 
  if(p1_inside && p2_inside) {
    return false;
  }
 
  if(p1_inside != p2_inside) {
    return true;
  }
 
  // at this point both end points of the line are outside the sphere
 
  // create a line and sphere bounding box
  MFloat bbox_line[2 * nDim];
  MFloat bbox_sphere[2 * nDim];
  for(MInt dim = 0; dim < nDim; dim++) {
    bbox_sphere[dim] = center[dim] - radius;
    bbox_sphere[dim + nDim] = center[dim] + radius;
    if(line[dim] < line[dim + nDim]) {
      bbox_line[dim] = line[dim];
      bbox_line[dim + nDim] = line[dim + nDim];
    } else {
      bbox_line[dim] = line[dim + nDim];
      bbox_line[dim + nDim] = line[dim];
    }
  }
 
  // check bounding boxes for overlap
  // a separation axis can be found in the cartesian axis
  if(!maia::geom::doBoxesOverlap<MFloat, nDim>(bbox_line, bbox_sphere)) {
    return false;
  }
 
  MFloat a = F0;
  MFloat b = F0;
  MFloat c = -POW2(radius);
  for(MInt dim = 0; dim < nDim; dim++) {
    const MFloat dl = line[dim + nDim] - line[dim];
    const MFloat dlc = line[dim] - center[dim];
    a += POW2(dl);
    b += dl * dlc;
    c += POW2(dlc);
  }
 
  const MFloat D = POW2(b / a) - c / a;
  if(D <= F0) {
    return false;
  }
 
  // at this point there have to be two valid intersections or none
  // it should be enough to check one of them
 
  const MFloat t_p = b / a + sqrt(D);
 
  // check the intersection if they lie on the line segment
  return (t_p <= F1) && (t_p >= F0);
}

isBoxInsideBox()

template<typename T , MInt nDim>

MBool maia::geom::isBoxInsideBox	(	const T *const	a,
		const T *const	b
	)

Definition at line 29 of file geometryanalytic.h.

                                                         {
  for(MInt dim = 0; dim < nDim; dim++) {
    if((a[dim] < b[dim]) || a[dim + nDim] > b[dim + nDim]) {
      return false;
    }
  }
  return true;
}

isBoxInsideSphere()

template<MInt nDim>

MBool maia::geom::isBoxInsideSphere	(	const MFloat *const	p,
		const MFloat *const	c,
		const MFloat	R
	)

Definition at line 96 of file geometryanalytic.h.

                                                                                      {
  MFloat pp[nDim];
  const MInt cornerCode[8][3] = {{0, 0, 0},    {nDim, 0, 0},    {0, nDim, 0},    {nDim, nDim, 0},
                                 {0, 0, nDim}, {nDim, 0, nDim}, {0, nDim, nDim}, {nDim, nDim, nDim}};
  for(MInt corner = 0; corner < IPOW2(nDim); corner++) {
    for(MInt dim = 0; dim < nDim; dim++) {
      pp[dim] = p[dim + cornerCode[corner][dim]];
    }
    if(!(maia::geom::isPointInsideSphere<nDim>(pp, c, R))) {
      return false;
    }
  }
  return true;
}

isPointInsideBox()

template<typename T , MInt nDim>

MBool maia::geom::isPointInsideBox	(	const T *const	p,
		const T *const	b
	)

Definition at line 39 of file geometryanalytic.h.

                                                           {
  for(MInt dim = 0; dim < nDim; dim++) {
    if((p[dim] < b[dim]) || p[dim] > b[dim + nDim]) {
      return false;
    }
  }
  return true;
}

isPointInsideRingSegment()

template<MInt nDim>

MBool maia::geom::isPointInsideRingSegment	(	const MFloat *const	p,
		const MFloat *const	c,
		const MFloat *const	R,
		const MFloat *const	phi,
		const MInt	axis
	)

Definition at line 62 of file geometryanalytic.h.

                                                                         {
  IF_CONSTEXPR(nDim == 2) {
    std::cerr << "H-Patch not implemented for 2D, not refining any cells! Take care with the angle calculation in 2D "
                 "when fixing this!"
              << std::endl;
    return false;
  }
 
  std::array<MFloat, nDim> normalVector;
  MFloat distance = F0;
 
  // compute normal vector from axis to current point
  for(MInt d = 0; d < nDim; d++) {
    if(d == axis) {
      normalVector[d] = F0;
    } else {
      normalVector[d] = p[d] - c[d];
      distance += normalVector[d] * normalVector[d];
    }
  }
 
  // check if cell is within r_min and r_max
  distance = sqrt(distance);
  if(distance < R[0] || distance > R[1]) return false;
 
  // check if cell is within phi_min and phi_max
  MFloat angle = atan2(normalVector[(axis + 2) % nDim], normalVector[(axis + 1) % nDim]) * 180.0 / PI;
  if(angle < phi[0] || angle > phi[1]) return false;
 
  return true;
}

isPointInsideSphere()

template<MInt nDim>

MBool maia::geom::isPointInsideSphere	(	const MFloat *const	p,
		const MFloat *const	c,
		const MFloat	R
	)

Definition at line 49 of file geometryanalytic.h.

                                                                                        {
  MFloat r = F0;
  for(MInt dim = 0; dim < nDim; dim++) {
    r += POW2(p[dim] - c[dim]);
  }
  return std::sqrt(r) < R;
}

isSphereInsideBox()

template<typename T , MInt nDim>

MBool maia::geom::isSphereInsideBox	(	const T *const	c,
		const T	R,
		const T *const	p
	)

Definition at line 112 of file geometryanalytic.h.

                                                                       {
  T bbox[2 * nDim];
  for(MInt dim = 0; dim < nDim; dim++) {
    bbox[dim] = c[dim] - R;
    bbox[dim + nDim] = c[dim] + R;
  }
  return maia::geom::isBoxInsideBox<T, nDim>(bbox, p);
}