110 logical :: initialized = .false.
111 real(kind=
dp) :: box_min(3) = huge(0.0_dp)
112 real(kind=
dp) :: box_max(3) = -huge(0.0_dp)
113 real(kind=
dp) :: center(3) = 0.0_dp
114 real(kind=
dp) :: diameter = huge(0.0_dp)
115 real(kind=
dp) :: surface_area = 0.0_dp
136 generic ::
operator(.lt.) => less
137 generic ::
operator(.gt.) => greater
184 class(*),
intent(in) :: object
185 real(kind=
dp),
intent(in),
optional :: padding
207 call neko_error(
"get_aabb: Unsupported object type")
210 if (
present(padding))
then
211 call box%add_padding(padding)
222 class(
aabb_t),
intent(inout) :: this
223 real(kind=
dp),
intent(in) :: padding
224 real(kind=
dp) :: box_min(3), box_max(3)
226 box_min = this%box_min - padding * (this%diameter)
227 box_max = this%box_max + padding * (this%diameter)
229 call this%init(box_min, box_max)
247 real(kind=
dp) :: box_min(3), box_max(3)
249 box_min = huge(0.0_dp); box_max = -huge(0.0_dp)
251 do i = 1, object%n_points()
253 box_min = min(box_min, pi%x)
254 box_max =
max(box_max, pi%x)
257 call box%init(box_min, box_max)
269 type(
mesh_t),
intent(in) :: object
275 do i = 1, object%nelv
276 temp_box =
get_aabb(object%elements(i))
277 box =
merge(box, temp_box)
297 do i = 1, object%nelv
299 box =
merge(box, temp_box)
319 do i = 1, object%nelv
321 box =
merge(box, temp_box)
334 subroutine aabb_init(this, lower_left_front, upper_right_back)
335 class(
aabb_t),
intent(inout) :: this
336 real(kind=
dp),
dimension(3),
intent(in) :: lower_left_front
337 real(kind=
dp),
dimension(3),
intent(in) :: upper_right_back
339 this%initialized = .true.
340 this%box_min = lower_left_front
341 this%box_max = upper_right_back
342 this%center = (this%box_min + this%box_max) / 2.0_dp
343 this%diameter = norm2(this%box_max - this%box_min)
345 this%surface_area = this%calculate_surface_area()
355 class(
aabb_t),
intent(in) :: this
356 real(kind=
dp),
dimension(3) :: min
363 class(
aabb_t),
intent(in) :: this
364 real(kind=
dp),
dimension(3) ::
max
371 class(
aabb_t),
intent(in) :: this
372 real(kind=
dp) :: width
374 width = this%box_max(1) - this%box_min(1)
379 class(
aabb_t),
intent(in) :: this
380 real(kind=
dp) :: depth
382 depth = this%box_max(2) - this%box_min(2)
387 class(
aabb_t),
intent(in) :: this
388 real(kind=
dp) :: height
390 height = this%box_max(3) - this%box_min(3)
395 class(
aabb_t),
intent(in) :: this
396 real(kind=
dp) :: diameter
398 diameter = this%diameter
403 class(
aabb_t),
intent(in) :: this
404 real(kind=
dp) :: surface_area
406 surface_area = this%surface_area
411 class(
aabb_t),
intent(in) :: this
412 real(kind=
dp),
dimension(3) :: center
419 class(
aabb_t),
intent(in) :: this
420 real(kind=
dp),
dimension(3) :: diagonal
422 diagonal = this%box_max - this%box_min
431 class(
aabb_t),
intent(in) :: this
432 class(
aabb_t),
intent(in) :: other
433 logical :: is_overlapping
435 if (.not. this%initialized .or. .not. other%initialized)
then
437 is_overlapping = .false.
440 is_overlapping = all(this%box_min .le. other%box_max) .and. &
441 all(this%box_max .ge. other%box_min)
448 class(
aabb_t),
intent(in) :: this
449 class(
aabb_t),
intent(in) :: other
450 logical :: is_contained
456 is_contained = all(this%box_min .le. other%box_min) .and. &
457 all(this%box_max .ge. other%box_max)
463 class(
aabb_t),
intent(in) :: this
464 real(kind=
dp),
dimension(3),
intent(in) :: p
465 logical :: is_contained
471 is_contained = all(p .ge. this%box_min) .and. all(p .le. this%box_max)
480 class(
aabb_t),
intent(in) :: box1
481 class(
aabb_t),
intent(in) :: box2
484 real(kind=
dp),
dimension(3) :: box_min, box_max
486 box_min = min(box1%box_min, box2%box_min)
487 box_max =
max(box1%box_max, box2%box_max)
489 call merged%init(box_min, box_max)
494 class(
aabb_t),
intent(in) :: box1
495 class(
aabb_t),
intent(in) :: box2
496 type(
aabb_t) :: intersected
498 real(kind=
dp),
dimension(3) :: box_min, box_max
500 box_min =
max(box1%box_min, box2%box_min)
501 box_max = min(box1%box_max, box2%box_max)
503 call intersected%init(box_min, box_max)
512 class(
aabb_t),
intent(in) :: this
513 real(kind=
dp) :: surface_area
515 surface_area = 2.0 * (&
516 & this%get_width() * this%get_height() &
517 & + this%get_width() * this%get_depth() &
518 & + this%get_height() * this%get_depth() &
528 class(
aabb_t),
intent(in) :: this
529 class(
aabb_t),
intent(in) :: other
533 if (.not. this%initialized .or. .not. other%initialized)
then
538 aabb_less = this%box_min(1) .lt. other%box_min(1)
539 equal = this%box_min(1) .le. other%box_min(1)
542 aabb_less = this%box_min(2) .lt. other%box_min(2)
543 equal = this%box_min(2) .le. other%box_min(2)
546 aabb_less = this%box_min(3) .lt. other%box_min(3)
547 equal = this%box_min(3) .le. other%box_min(3)
550 aabb_less = this%box_max(1) .lt. other%box_max(1)
551 equal = this%box_max(1) .le. other%box_max(1)
554 aabb_less = this%box_max(2) .lt. other%box_max(2)
555 equal = this%box_max(2) .le. other%box_max(2)
558 aabb_less = this%box_max(3) .lt. other%box_max(3)
565 class(
aabb_t),
intent(in) :: this
566 class(
aabb_t),
intent(in) :: other
570 if (.not. this%initialized .or. .not. other%initialized)
then
576 equal = this%box_min(1) .ge. other%box_min(1)
580 equal = this%box_min(2) .ge. other%box_min(2)
584 equal = this%box_min(3) .ge. other%box_min(3)
588 equal = this%box_max(1) .ge. other%box_max(1)
592 equal = this%box_max(2) .ge. other%box_max(2)
Axis Aligned Bounding Box (aabb) implementation in Fortran.
pure real(kind=dp) function, dimension(3) aabb_get_max(this)
Get the maximum point of the aabb.
pure real(kind=dp) function aabb_get_surface_area(this)
Get the surface area of the aabb.
pure real(kind=dp) function aabb_get_diameter(this)
Get the diameter length of the aabb.
subroutine add_padding(this, padding)
Add padding to the aabb.
logical function aabb_contains_point(this, p)
Check if this aabb contains a point.
subroutine aabb_init(this, lower_left_front, upper_right_back)
Initialize the aabb.
type(aabb_t) function get_aabb_tri_mesh(object)
Get the aabb of a triangular mesh.
logical function aabb_overlaps(this, other)
Check if two aabbs are overlapping.
type(aabb_t) function get_aabb_mesh(object)
Get the aabb of a mesh.
type(aabb_t) function intersection_aabb(box1, box2)
Get the intersection of two aabbs.
pure real(kind=dp) function aabb_get_height(this)
Get the height of the aabb. Also known as the z-axis length.
pure logical function aabb_greater(this, other)
Greater than comparison operator.
type(aabb_t) function, public get_aabb(object, padding)
Construct the aabb of a predefined object.
logical function aabb_contains_other(this, other)
Check if this aabb contains another aabb.
pure real(kind=dp) function calculate_surface_area(this)
Calculate the surface area of the aabb.
pure logical function aabb_less(this, other)
Less than comparison operator.
type(aabb_t) function get_aabb_element(object)
Get the aabb of an arbitrary element.
type(aabb_t) function get_aabb_tet_mesh(object)
Get the aabb of a tetrahedral mesh.
pure real(kind=dp) function, dimension(3) aabb_get_center(this)
Get the center of the aabb.
pure real(kind=dp) function aabb_get_width(this)
Get the width of the aabb. Also known as the x-axis length.
pure real(kind=dp) function, dimension(3) aabb_get_diagonal(this)
Get the diagonal of the aabb.
type(aabb_t) function merge_aabb(box1, box2)
Merge two aabbs.
pure real(kind=dp) function aabb_get_depth(this)
Get the depth of the aabb. Also known as the y-axis length.
pure real(kind=dp) function, dimension(3) aabb_get_min(this)
Get the minimum point of the aabb.
Defines a hexahedron element.
integer, parameter, public dp
Defines a quadrilateral element.
Defines a tetrahedral mesh.
Defines a tetrahedral element.
Defines a triangular surface mesh.
Defines a triangular element.
Axis Aligned Bounding Box (aabb) data structure.
Base type for an element.
A point in with coordinates .