Skip to content

Commit

Permalink
Add generic function for compute multiplicity.
Browse files Browse the repository at this point in the history 
- Add compute_multiplicity2.
- Update and fix NURBS functions.
- Remove unused variables.
  • Loading branch information
@gha3mi
gha3mi committed Apr 6, 2024
1 parent 90f2e25 commit 35de96f
Show file tree
Hide file tree
Showing 4 changed files with 68 additions and 46 deletions.
11 changes: 4 additions & 7 deletions src/NURBS/forcad_nurbs_curve.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -496,15 +496,13 @@ pure subroutine insert_knots(this,Xth,r)
integer, intent(in) :: r(:)
integer :: k, i, s, dim, j, n_new
real(rk), allocatable :: Xcw(:,:), Xcw_new(:,:), Xc_new(:,:), Wc_new(:), knot_new(:)
integer, allocatable :: mlp(:)

if (allocated(this%Wc)) then ! NURBS

do i = 1, size(Xth)
k = findspan(this%nc-1,this%order,Xth(i),this%knot)
if (this%knot(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot)
s = mlp(k)
if (this%knot(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot,Xth(i))
else
s = 0
end if
Expand Down Expand Up @@ -543,9 +541,8 @@ pure subroutine insert_knots(this,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc-1,this%order,Xth(i),this%knot)
if (this%knot(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot)
s = mlp(k)
if (this%knot(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot,Xth(i))
else
s = 0
end if
Expand Down
23 changes: 9 additions & 14 deletions src/NURBS/forcad_nurbs_surface.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -722,9 +722,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)
integer, intent(in) :: dir
real(rk), intent(in) :: Xth(:)
integer, intent(in) :: r(:)
integer :: k, i, s, dim, j, n_new, jj
integer :: k, i, s, dim, j, n_new
real(rk), allocatable :: Xc(:,:), Xcw(:,:), Xcw_new(:,:), Xc_new(:,:), Wc_new(:), knot_new(:)
integer, allocatable :: mlp(:)
real(rk), allocatable:: Xc3(:,:,:)


Expand All @@ -734,9 +733,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc(1)-1,this%order(1),Xth(i),this%knot1)
if (this%knot1(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot1)
s = mlp(k)
if (this%knot1(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot1, Xth(i))
else
s = 0
end if
Expand Down Expand Up @@ -780,9 +778,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc(1)-1,this%order(1),Xth(i),this%knot1)
if (this%knot1(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot1)
s = mlp(k)
if (this%knot1(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot1, Xth(i))
else
s = 0
end if
Expand Down Expand Up @@ -820,9 +817,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc(2)-1,this%order(2),Xth(i),this%knot2)
if (this%knot2(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot2)
s = mlp(k)
if (this%knot2(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot2, Xth(i))
else
s = 0
end if
Expand Down Expand Up @@ -871,9 +867,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc(2)-1,this%order(2),Xth(i),this%knot2)
if (this%knot2(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot2)
s = mlp(k)
if (this%knot2(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot2, Xth(i))
else
s = 0
end if
Expand Down
38 changes: 15 additions & 23 deletions src/NURBS/forcad_nurbs_volume.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -816,10 +816,9 @@ pure subroutine insert_knots(this, dir ,Xth,r)
integer, intent(in) :: dir
real(rk), intent(in) :: Xth(:)
integer, intent(in) :: r(:)
integer :: k, i, s, dim, j, n_new, jj
integer :: k, i, s, dim, j, n_new
real(rk), allocatable :: Xc(:,:), Xcw(:,:), Xcw_new(:,:), Xc_new(:,:), Wc_new(:), knot_new(:)
integer, allocatable :: mlp(:)
real(rk), allocatable:: Xc3(:,:,:), Xc4(:,:,:,:)
real(rk), allocatable :: Xc4(:,:,:,:)


if (dir == 1) then ! direction 1
Expand All @@ -828,9 +827,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc(1)-1,this%order(1),Xth(i),this%knot1)
if (this%knot1(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot1)
s = mlp(k)
if (this%knot1(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot1, Xth(i))
else
s = 0
end if
Expand Down Expand Up @@ -874,9 +872,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc(1)-1,this%order(1),Xth(i),this%knot1)
if (this%knot1(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot1)
s = mlp(k)
if (this%knot1(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot1, Xth(i))
else
s = 0
end if
Expand Down Expand Up @@ -914,9 +911,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc(2)-1,this%order(2),Xth(i),this%knot2)
if (this%knot2(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot2)
s = mlp(k)
if (this%knot2(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot2, Xth(i))
else
s = 0
end if
Expand Down Expand Up @@ -965,9 +961,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc(2)-1,this%order(2),Xth(i),this%knot2)
if (this%knot2(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot2)
s = mlp(k)
if (this%knot2(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot2, Xth(i))
else
s = 0
end if
Expand Down Expand Up @@ -1009,9 +1004,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc(3)-1,this%order(3),Xth(i),this%knot3)
if (this%knot3(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot3)
s = mlp(k)
if (this%knot3(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot3, Xth(i))
else
s = 0
end if
Expand Down Expand Up @@ -1060,9 +1054,8 @@ pure subroutine insert_knots(this, dir ,Xth,r)

do i = 1, size(Xth)
k = findspan(this%nc(3)-1,this%order(3),Xth(i),this%knot3)
if (this%knot3(k) == Xth(i)) then
mlp = compute_multiplicity(this%knot3)
s = mlp(k)
if (this%knot3(k+1) == Xth(i)) then
s = compute_multiplicity(this%knot3, Xth(i))
else
s = 0
end if
Expand Down Expand Up @@ -1115,8 +1108,7 @@ pure subroutine elevate_degree(this, dir, t)
integer, intent(in) :: t
real(rk), allocatable :: Xc(:,:), Xcw(:,:), Xcw_new(:,:), Xc_new(:,:), Wc_new(:), knot_new(:)
integer :: nc_new, dim, j
integer, allocatable :: mlp(:)
real(rk), allocatable:: Xc3(:,:,:), Xc4(:,:,:,:)
real(rk), allocatable:: Xc4(:,:,:,:)


if (dir == 1) then ! direction 1
Expand Down
42 changes: 40 additions & 2 deletions src/forcad_utils.f90
View file
Original file line number Diff line number Diff line change
Expand Up @@ -24,6 +24,14 @@ module forcad_utils
end interface
!===============================================================================


!===============================================================================
interface compute_multiplicity
module procedure compute_multiplicity1
module procedure compute_multiplicity2
end interface
!===============================================================================

contains

!===============================================================================
Expand Down Expand Up @@ -69,7 +77,7 @@ pure function basis_bspline_der(Xt, knot, nc, order) result(dB)
real(rk), allocatable :: dB(:)
real(rk), allocatable :: Nt(:,:), dNt_dXt(:,:)
real(rk) :: R, L, Rp, Lp, knot_i, knot_ip, knot_jk, knot_jkm, knot_end, a, b, c, d
integer :: i, p, k, n, m, jk
integer :: i, k, n, m, jk

k = order + 1
n = nc - 1
Expand Down Expand Up @@ -304,7 +312,7 @@ pure function cmp_elemConn_C0_V(nnode1,nnode2,nnode3,p1,p2,p3) result(elemConn)
!===============================================================================
!> author: Seyed Ali Ghasemi
!> license: BSD 3-Clause
pure function compute_multiplicity(knot) result(multiplicity)
pure function compute_multiplicity1(knot) result(multiplicity)
real(rk), intent(in) :: knot(:)
integer, dimension(:), allocatable :: multiplicity
integer :: i, count
Expand All @@ -331,6 +339,36 @@ pure function compute_multiplicity(knot) result(multiplicity)
!===============================================================================


!===============================================================================
!> author: Seyed Ali Ghasemi
!> license: BSD 3-Clause
pure function compute_multiplicity2(knot, Xth) result(multiplicity)
real(rk), dimension(:), intent(in) :: knot
real(rk), intent(in) :: Xth
integer :: multiplicity
integer :: i, count, size_knot

size_knot = size(knot)
multiplicity = 0
i = 1
do while (i <= size_knot)
if (knot(i) == Xth) then
count = 1
do while (i + count <= size_knot .and. knot(i + count) == Xth)
count = count + 1
end do
if (count > multiplicity) then
multiplicity = count
end if
i = i + count
else
i = i + 1
end if
end do
end function
!===============================================================================


!===============================================================================
!> author: Seyed Ali Ghasemi
!> license: BSD 3-Clause
Expand Down

0 comments on commit 35de96f

Please sign in to comment.