Homec4science
Diffusion SPClibs (feature/refactoring) Browse

spclibs/pppack/knots.f90
3324689e437cfeature/refactoring

Last Change
Raw File

knots.f90
View Options

!>
!> @file knots.f90
!>
!> @brief
!>
!> @copyright
!> Copyright (©) 2021 EPFL (Ecole Polytechnique Fédérale de Lausanne)
!> SPC (Swiss Plasma Center)
!>
!> SPClibs is free software: you can redistribute it and/or modify it under
!> the terms of the GNU Lesser General Public License as published by the Free
!> Software Foundation, either version 3 of the License, or (at your option)
!> any later version.
!>
!> SPClibs is distributed in the hope that it will be useful, but WITHOUT ANY
!> WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
!> FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
!>
!> You should have received a copy of the GNU Lesser General Public License
!> along with this program. If not, see <https://www.gnu.org/licenses/>.
!>
!> @author
!> (in alphabetical order)
!> @author Trach-Minh Tran <trach-minh.tran@epfl.ch>
!>
subroutine knots ( break, l, kpm, t, n )
!*************************************************************************
!
!! KNOTS is to be called in COLLOC.
!
! Discussion:
!
! From the given breakpoint sequence BREAK the routine constructs the
! knot sequence T so that
!
! SPLINE(K+M,T) = PP(K+M,BREAK)
!
! with M-1 continuous derivatives. This means that
!
! t(1),...,t(n+kpm) = break(1) kpm times, then break(2),...,
! break(l) each k times, then, finally, break(l+1) kpm times.
!
! Reference:
!
! Carl DeBoor,
! A Practical Guide to Splines,
! Springer Verlag.
!
! Parameters:
!
! Input, real ( kind = 8 ) BREAK(L+1), the breakpoint sequence.
!
! Input, integer L, the number of intervals or pieces.
!
! Input, integer KPM, = K+M, the order of the piecewise polynomial
! function or spline.
!
! Output, real ( kind = 8 ) T(N+KPM), the knot sequence.
!
! Output, integer N, = L*K+M = the dimension of SPLINE(K+M,T).
!
implicit none
integer kpm
integer l
integer n
real ( kind = 8 ) break(l+1)
integer iside
integer j
integer jj
integer jjj
integer k
integer ll
integer m
real ( kind = 8 ) t(*)
real ( kind = 8 ) xside
common /side/ m,iside,xside(10)
k = kpm-m
n = l*k+m
jj = n+kpm
jjj = l+1
do ll = 1, kpm
t(jj) = break(jjj)
jj = jj-1
end do
do j = 1, l
jjj = jjj-1
do ll = 1, k
t(jj) = break(jjj)
jj = jj-1
end do
end do
t(1:kpm) = break(1)
return
end
c4science · Help