55!
66! # Description
77!
8- ! Multidimensional (2D -6D) B-Spline interpolation of data on a regular grid.
8+ ! Multidimensional (1D -6D) B-Spline interpolation of data on a regular grid.
99! Basic subroutine interface.
1010!
1111! # Notes
1414! The original Fortran 77 routines were converted to free-form source.
1515! Some of them are relatively unchanged from the originals, but some have
1616! been extensively refactored. In addition, new routines for
17- ! 4d, 5d, and 6d interpolation were also created (these are simply
17+ ! 1d, 4d, 5d, and 6d interpolation were also created (these are simply
1818! extensions of the same algorithm into higher dimensions).
1919!
2020! # See also
@@ -48,6 +48,7 @@ module bspline_sub_module
4848 integer ,parameter :: wp = real64 ! ! Real precision
4949
5050 ! main routines:
51+ public :: db1ink, db1val
5152 public :: db2ink, db2val
5253 public :: db3ink, db3val
5354 public :: db4ink, db4val
@@ -58,6 +59,134 @@ module bspline_sub_module
5859
5960! *****************************************************************************************
6061
62+ ! *****************************************************************************************
63+ ! > Determines the parameters of a function that interpolates
64+ ! the one-dimensional gridded data
65+ ! $$ [x(i),\mathrm{fcn}(i)] ~\mathrm{for}~ i=1,..,n_x $$
66+ ! The interpolating function and its derivatives may
67+ ! subsequently be evaluated by the function [[db1val]].
68+ !
69+ ! # History
70+ !
71+ ! * Jacob Williams, 10/30/2015 : Created 1D routine.
72+
73+ subroutine db1ink (x ,nx ,fcn ,kx ,tx ,bcoef ,iflag )
74+
75+ implicit none
76+
77+ integer ,intent (in ) :: nx ! ! Number of x abcissae
78+ integer ,intent (in ) :: kx ! ! The order of spline pieces in x (>= 2, < nx). (order = polynomial degree + 1)
79+ real (wp),dimension (nx),intent (in ) :: x ! ! Array of x abcissae. Must be strictly increasing.
80+ real (wp),dimension (nx),intent (in ) :: fcn ! ! Array of function values to interpolate. fcn(i) should
81+ ! ! contain the function value at the point x(i)
82+ real (wp),dimension (nx+ kx),intent (inout ) :: tx ! ! The knots in the x direction for the spline interpolant.
83+ ! ! If iflag=0 these are chosen by db1ink.
84+ ! ! If iflag=1 these are specified by the user.
85+ ! ! Must be non-decreasing.
86+ real (wp),dimension (nx),intent (out ) :: bcoef ! ! Array of coefficients of the b-spline interpolant.
87+ integer ,intent (inout ) :: iflag ! ! **on input:** 0 = knot sequence chosen by db1ink.
88+ ! ! 1 = knot sequence chosen by user.
89+ ! ! **on output:** 1 = successful execution.
90+ ! ! 2 = iflag out of range.
91+ ! ! 3 = nx out of range.
92+ ! ! 4 = kx out of range.
93+ ! ! 5 = x not strictly increasing.
94+ ! ! 6 = tx not non-decreasing.
95+
96+
97+ real (wp),dimension (nx) :: temp
98+ real (wp),dimension (2 * kx* (nx+1 )) :: work
99+ logical :: status_ok
100+
101+ ! check validity of inputs
102+
103+ call check_inputs(' db1ink'
104+ iflag,&
105+ nx= nx,&
106+ kx= kx,&
107+ x= x,&
108+ tx= tx,&
109+ status_ok= status_ok)
110+
111+ if (status_ok) then
112+
113+ ! choose knots
114+
115+ if (iflag == 0 ) then
116+ call dbknot(x,nx,kx,tx)
117+ end if
118+
119+ ! construct b-spline coefficients
120+
121+ call dbtpcf(x,nx,fcn,nx,1 ,tx,kx,bcoef,work)
122+
123+ iflag = 1
124+
125+ end if
126+
127+ end subroutine db1ink
128+ ! *****************************************************************************************
129+
130+ ! *****************************************************************************************
131+ ! > Evaluates the tensor product piecewise polynomial
132+ ! interpolant constructed by the routine [[db1ink]] or one of its
133+ ! derivatives at the point xval.
134+ !
135+ ! To evaluate the interpolant itself, set idx=0,
136+ ! to evaluate the first partial with respect to x, set idx=1, and so on.
137+ !
138+ ! db1val returns 0.0 if (xval,yval) is out of range. that is, if
139+ ! ```fortran
140+ ! xval < tx(1) .or. xval > tx(nx+kx)
141+ ! ```
142+ ! if the knots tx were chosen by [[db1ink]], then this is equivalent to:
143+ ! ```fortran
144+ ! xval < x(1) .or. xval > x(nx)+epsx
145+ ! ```
146+ ! where
147+ ! ```fortran
148+ ! epsx = 0.1*(x(nx)-x(nx-1))
149+ ! ```
150+ !
151+ ! The input quantities tx, nx, kx, and bcoef should be
152+ ! unchanged since the last call of [[db1ink]].
153+ !
154+ ! # History
155+ !
156+ ! * Jacob Williams, 10/30/2015 : Created 1D routine.
157+
158+ subroutine db1val (xval ,idx ,tx ,nx ,kx ,bcoef ,f ,iflag ,inbvx_in )
159+
160+ implicit none
161+
162+ integer ,intent (in ) :: idx ! ! x derivative of piecewise polynomial to evaluate.
163+ integer ,intent (in ) :: nx ! ! the number of interpolation points in x. (same as in last call to db1ink)
164+ integer ,intent (in ) :: kx ! ! order of polynomial pieces in x. (same as in last call to db1ink)
165+ real (wp),intent (in ) :: xval ! ! x coordinate of evaluation point.
166+ real (wp),dimension (nx+ kx),intent (in ) :: tx ! ! sequence of knots defining the piecewise polynomial in the x direction. (same as in last call to db1ink)
167+ real (wp),dimension (nx),intent (in ) :: bcoef ! ! the b-spline coefficients computed by db1ink.
168+ real (wp),intent (out ) :: f ! ! interpolated value
169+ integer ,intent (out ) :: iflag ! ! status flag: 0 : no errors, /=0 : error
170+ integer ,intent (in ),optional :: inbvx_in ! ! inbvx initialization parameter from class
171+
172+ real (wp),dimension (3 * kx) :: work
173+
174+ integer ,save :: inbvx = 1
175+
176+ if (present (inbvx_in)) inbvx = inbvx_in
177+
178+ f = 0.0_wp
179+ iflag = 1
180+
181+ if (xval< tx(1 ) .or. xval> tx(nx+ kx)) return
182+
183+ iflag = 0
184+
185+ f = dbvalu(tx,bcoef,nx,kx,idx,xval,inbvx,work)
186+
187+ end subroutine db1val
188+ ! *****************************************************************************************
189+
61190! *****************************************************************************************
62191! > Determines the parameters of a function that interpolates
63192! the two-dimensional gridded data
0 commit comments