Skip to content

Commit

Permalink
Adding support for degree=4 case for Fortran surface subdivide_nodes.
Browse files Browse the repository at this point in the history
  • Loading branch information
dhermes committed Nov 2, 2017
1 parent 8beb1ac commit 0b2b1f3
Show file tree
Hide file tree
Showing 2 changed files with 131 additions and 2 deletions.
124 changes: 124 additions & 0 deletions src/bezier/surface.f90
Expand Up @@ -391,6 +391,130 @@ subroutine subdivide_nodes( &
nodes_d(8, :) = 0.5_dp * (nodes(8, :) + nodes(10, :))
nodes_d(9, :) = 0.5_dp * (nodes(9, :) + nodes(10, :))
nodes_d(10, :) = nodes(10, :)
else if (degree == 4) then
nodes_a(1, :) = nodes(1, :)
nodes_a(2, :) = 0.5_dp * (nodes(1, :) + nodes(2, :))
nodes_a(3, :) = 0.25_dp * (nodes(1, :) + 2 * nodes(2, :) + nodes(3, :))
nodes_a(4, :) = 0.125_dp * ( &
nodes(1, :) + 3 * nodes(2, :) + 3 * nodes(3, :) + nodes(4, :))
nodes_a(5, :) = 0.0625_dp * ( &
nodes(1, :) + 4 * nodes(2, :) + 6 * nodes(3, :) + &
4 * nodes(4, :) + nodes(5, :))
nodes_a(6, :) = 0.5_dp * (nodes(1, :) + nodes(6, :))
nodes_a(7, :) = 0.25_dp * ( &
nodes(1, :) + nodes(2, :) + nodes(6, :) + nodes(7, :))
nodes_a(8, :) = 0.125_dp * ( &
nodes(1, :) + 2 * nodes(2, :) + nodes(3, :) + nodes(6, :) + &
2 * nodes(7, :) + nodes(8, :))
nodes_a(9, :) = 0.0625_dp * ( &
nodes(1, :) + 3 * nodes(2, :) + 3 * nodes(3, :) + nodes(4, :) + &
nodes(6, :) + 3 * nodes(7, :) + 3 * nodes(8, :) + nodes(9, :))
nodes_a(10, :) = 0.25_dp * ( &
nodes(1, :) + 2 * nodes(6, :) + nodes(10, :))
nodes_a(11, :) = 0.125_dp * ( &
nodes(1, :) + nodes(2, :) + 2 * nodes(6, :) + 2 * nodes(7, :) + &
nodes(10, :) + nodes(11, :))
nodes_a(12, :) = 0.0625_dp * ( &
nodes(1, :) + 2 * nodes(2, :) + nodes(3, :) + 2 * nodes(6, :) + &
4 * nodes(7, :) + 2 * nodes(8, :) + nodes(10, :) + &
2 * nodes(11, :) + nodes(12, :))
nodes_a(13, :) = 0.125_dp * ( &
nodes(1, :) + 3 * nodes(6, :) + 3 * nodes(10, :) + nodes(13, :))
nodes_a(14, :) = 0.0625_dp * ( &
nodes(1, :) + nodes(2, :) + 3 * nodes(6, :) + 3 * nodes(7, :) + &
3 * nodes(10, :) + 3 * nodes(11, :) + nodes(13, :) + nodes(14, :))
nodes_a(15, :) = 0.0625_dp * ( &
nodes(1, :) + 4 * nodes(6, :) + 6 * nodes(10, :) + &
4 * nodes(13, :) + nodes(15, :))
nodes_b(1, :) = 0.0625_dp * ( &
nodes(5, :) + 4 * nodes(9, :) + 6 * nodes(12, :) + &
4 * nodes(14, :) + nodes(15, :))
nodes_b(2, :) = 0.0625_dp * ( &
nodes(4, :) + 3 * nodes(8, :) + nodes(9, :) + 3 * nodes(11, :) + &
3 * nodes(12, :) + nodes(13, :) + 3 * nodes(14, :) + nodes(15, :))
nodes_b(3, :) = 0.0625_dp * ( &
nodes(3, :) + 2 * nodes(7, :) + 2 * nodes(8, :) + nodes(10, :) + &
4 * nodes(11, :) + nodes(12, :) + 2 * nodes(13, :) + &
2 * nodes(14, :) + nodes(15, :))
nodes_b(4, :) = 0.0625_dp * ( &
nodes(2, :) + nodes(6, :) + 3 * nodes(7, :) + 3 * nodes(10, :) + &
3 * nodes(11, :) + 3 * nodes(13, :) + nodes(14, :) + nodes(15, :))
nodes_b(5, :) = nodes_a(15, :)
nodes_b(6, :) = 0.0625_dp * ( &
nodes(4, :) + nodes(5, :) + 3 * nodes(8, :) + 3 * nodes(9, :) + &
3 * nodes(11, :) + 3 * nodes(12, :) + nodes(13, :) + nodes(14, :))
nodes_b(7, :) = 0.0625_dp * ( &
nodes(3, :) + nodes(4, :) + 2 * nodes(7, :) + 3 * nodes(8, :) + &
nodes(9, :) + nodes(10, :) + 3 * nodes(11, :) + &
2 * nodes(12, :) + nodes(13, :) + nodes(14, :))
nodes_b(8, :) = 0.0625_dp * ( &
nodes(2, :) + nodes(3, :) + nodes(6, :) + 3 * nodes(7, :) + &
2 * nodes(8, :) + 2 * nodes(10, :) + 3 * nodes(11, :) + &
nodes(12, :) + nodes(13, :) + nodes(14, :))
nodes_b(9, :) = nodes_a(14, :)
nodes_b(10, :) = 0.0625_dp * ( &
nodes(3, :) + 2 * nodes(4, :) + nodes(5, :) + 2 * nodes(7, :) + &
4 * nodes(8, :) + 2 * nodes(9, :) + nodes(10, :) + &
2 * nodes(11, :) + nodes(12, :))
nodes_b(11, :) = 0.0625_dp * ( &
nodes(2, :) + 2 * nodes(3, :) + nodes(4, :) + nodes(6, :) + &
3 * nodes(7, :) + 3 * nodes(8, :) + nodes(9, :) + nodes(10, :) + &
2 * nodes(11, :) + nodes(12, :))
nodes_b(12, :) = nodes_a(12, :)
nodes_b(13, :) = 0.0625_dp * ( &
nodes(2, :) + 3 * nodes(3, :) + 3 * nodes(4, :) + nodes(5, :) + &
nodes(6, :) + 3 * nodes(7, :) + 3 * nodes(8, :) + nodes(9, :))
nodes_b(14, :) = nodes_a(9, :)
nodes_b(15, :) = nodes_a(5, :)
nodes_c(1, :) = nodes_a(5, :)
nodes_c(2, :) = 0.125_dp * ( &
nodes(2, :) + 3 * nodes(3, :) + 3 * nodes(4, :) + nodes(5, :))
nodes_c(3, :) = 0.25_dp * ( &
nodes(3, :) + 2 * nodes(4, :) + nodes(5, :))
nodes_c(4, :) = 0.5_dp * ( &
nodes(4, :) + nodes(5, :))
nodes_c(5, :) = nodes(5, :)
nodes_c(6, :) = nodes_b(13, :)
nodes_c(7, :) = 0.125_dp * ( &
nodes(3, :) + 2 * nodes(4, :) + nodes(5, :) + nodes(7, :) + &
2 * nodes(8, :) + nodes(9, :))
nodes_c(8, :) = 0.25_dp * ( &
nodes(4, :) + nodes(5, :) + nodes(8, :) + nodes(9, :))
nodes_c(9, :) = 0.5_dp * (nodes(5, :) + nodes(9, :))
nodes_c(10, :) = nodes_b(10, :)
nodes_c(11, :) = 0.125_dp * ( &
nodes(4, :) + nodes(5, :) + 2 * nodes(8, :) + 2 * nodes(9, :) + &
nodes(11, :) + nodes(12, :))
nodes_c(12, :) = 0.25_dp * ( &
nodes(5, :) + 2 * nodes(9, :) + nodes(12, :))
nodes_c(13, :) = nodes_b(6, :)
nodes_c(14, :) = 0.125_dp * ( &
nodes(5, :) + 3 * nodes(9, :) + 3 * nodes(12, :) + nodes(14, :))
nodes_c(15, :) = nodes_b(1, :)
nodes_d(1, :) = nodes_a(15, :)
nodes_d(2, :) = nodes_b(4, :)
nodes_d(3, :) = nodes_b(3, :)
nodes_d(4, :) = nodes_b(2, :)
nodes_d(5, :) = nodes_b(1, :)
nodes_d(6, :) = 0.125_dp * ( &
nodes(6, :) + 3 * nodes(10, :) + 3 * nodes(13, :) + nodes(15, :))
nodes_d(7, :) = 0.125_dp * ( &
nodes(7, :) + nodes(10, :) + 2 * nodes(11, :) + &
2 * nodes(13, :) + nodes(14, :) + nodes(15, :))
nodes_d(8, :) = 0.125_dp * ( &
nodes(8, :) + 2 * nodes(11, :) + nodes(12, :) + nodes(13, :) + &
2 * nodes(14, :) + nodes(15, :))
nodes_d(9, :) = 0.125_dp * ( &
nodes(9, :) + 3 * nodes(12, :) + 3 * nodes(14, :) + nodes(15, :))
nodes_d(10, :) = 0.25_dp * ( &
nodes(10, :) + 2 * nodes(13, :) + nodes(15, :))
nodes_d(11, :) = 0.25_dp * ( &
nodes(11, :) + nodes(13, :) + nodes(14, :) + nodes(15, :))
nodes_d(12, :) = 0.25_dp * ( &
nodes(12, :) + 2 * nodes(14, :) + nodes(15, :))
nodes_d(13, :) = 0.5_dp * (nodes(13, :) + nodes(15, :))
nodes_d(14, :) = 0.5_dp * (nodes(14, :) + nodes(15, :))
nodes_d(15, :) = nodes(15, :)
end if

end subroutine subdivide_nodes
Expand Down
9 changes: 7 additions & 2 deletions tests/fortran/test_surface.f90
Expand Up @@ -618,7 +618,7 @@ subroutine test_subdivide_nodes(success)
all(nodes_d2 == expected_d2))
call print_status(name, case_id, case_success, success)

! CASE 4: Evaluate subdivided parts of a linear surface.
! CASE 4: Evaluate subdivided parts of a quadratic surface.
call subdivide_points_check( &
6, 2, 2, 219803, 7086, case_success)
call print_status(name, case_id, case_success, success)
Expand Down Expand Up @@ -684,11 +684,16 @@ subroutine test_subdivide_nodes(success)
all(nodes_d3 == expected_d3))
call print_status(name, case_id, case_success, success)

! CASE 6: Evaluate subdivided parts of a linear surface.
! CASE 6: Evaluate subdivided parts of a cubic surface.
call subdivide_points_check( &
10, 2, 3, 439028340, 2184938, case_success)
call print_status(name, case_id, case_success, success)

! CASE 8: Evaluate subdivided parts of a quartic surface.
call subdivide_points_check( &
15, 2, 4, 1029038, 890012, case_success)
call print_status(name, case_id, case_success, success)

end subroutine test_subdivide_nodes

subroutine subdivide_points_check( &
Expand Down

0 comments on commit 0b2b1f3

Please sign in to comment.