DersOneBasisFun Subroutine
private pure subroutine DersOneBasisFun(me, i, u, n, ders)
Computes a single basis function and its derivatives based on Algorithm A2.5 in the NURBS Book [^1].
Implementation Details of Algorithm A2.5
By using Eq. (2.9), the algorithm computes for fixed where ; . The th derivative is returned in real :: ders(k). For example, if and , then
Using triangular tables, we must compute
The algorithm computes and stores the entire triangular table corresponding to ; to get the th derivative, loads the column of the table which contains the functions of degree , and compute the remaining portion of the triangle.
Note: This algorithm does not apply to rational basis functions.
Arguments
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| class(basis), | intent(in) | :: | me | |||
| integer, | intent(in) | :: | i | Basis function index |
||
| real(kind=wp), | intent(in) | :: | u | Given value |
||
| integer, | intent(in) | :: | n | Number of derivatives () |
||
| real(kind=wp), | intent(out) | :: | ders(0:n) |
|
Calls
Source Code
pure subroutine DersOneBasisFun(me,i,u,n,ders)
class(basis), intent(in) :: me
integer, intent(in) :: i !< Basis function index
real(wp), intent(in) :: u !< Given \( u \) value
integer, intent(in) :: n !< Number of derivatives (\(n \leq p\))
real(wp), intent(out) :: ders(0:n) !< \( = \{ N_{i,p}^{(0)}(u), \ldots, N_{i,p}^{(k)}(u) \} \)
!locals
integer :: m
integer :: j, jj, k
real(wp) :: ND(0:me%p), Ni(0:me%p,0:n), Uleft, Uright, saved, temp
associate( p => me%p, Ui => me%kv ) !get rid of % signs
! set highest knot vector index
m = me%nk-1
ND = 0.0_wp
Ni = 0.0_wp
ders = 0.0_wp
! if knot is outside of span range
if (u < Ui(i) .or. u >= Ui(i+p+1))then ! /* Local property */
do k=0,n
ders(k) = 0.0_wp
end do
return
end if
do j=0,p ! /* Initialize zeroth-degree functions */
if (u >= Ui(i+j) .and. u < Ui(i+j+1)) then
Ni(j,0) = 1.0_wp
else
Ni(j,0) = 0.0_wp
end if
end do
do k=1,p ! / * Compute full triangular table * /
! Detecting zeros saves computations
if (Ni(0,k-1) == 0.0_wp)then
saved = 0.0_wp
else
saved = ((u-Ui(i)) * Ni(0,k-1)) / (Ui(i+k)-Ui(i))
end if
do j=0,p-k
Uleft = Ui(i+j+1)
Uright = Ui(i+j+k+1)
! Zero detection
if (Ni(j+1,k-1) == 0.0_wp)then
Ni(j,k) = saved
saved = 0.0_wp
else
temp = Ni(j+1,k-1)/(Uright-Uleft)
Ni(j,k) = saved+(Uright-u)*temp
saved = (u-Uleft)*temp
end if
end do
end do
ders(0) = Ni(0,p) ! /* The function value */
do k=1,n ! /* Compute the derivatives */
ND = 0.0_wp
do j=0,k ! /* Load appropriate column */
ND(j) = Ni(j,p-k)
end do
do jj=1,k ! /* Compute table of width k */
if (ND(0) == 0.0_wp)then
saved = 0.0_wp
else
saved = ND(0) / (Ui(i+p-k+jj) - Ui(i))
end if
do j=0,k-jj
Uleft = Ui(i+j+1)
Uright = Ui(i+j+p-k+jj+1) !Wrong in The NURBS Book: -k is missing.
if (ND(j+1) == 0.0_wp)then
ND(j) = (p-k+jj)*saved
saved = 0.0_wp
else
temp = ND(j+1)/(Uright-Uleft);
ND(j) = (p-k+jj)*(saved-temp);
saved = temp;
end if
end do
end do
ders(k) = ND(0) ! /* kth derivative */
end do
end associate
end subroutine DersOneBasisFun