Determines the parameters of a function that interpolates the five-dimensional gridded data:
for:
The interpolating function and its derivatives may subsequently be evaluated by the function db5val.
See db3ink header for more details.
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
real(kind=wp), | intent(in), | dimension(:) | :: | x |
|
|
integer(kind=ip), | intent(in) | :: | nx |
number of abcissae ( ) |
||
real(kind=wp), | intent(in), | dimension(:) | :: | y |
|
|
integer(kind=ip), | intent(in) | :: | ny |
number of abcissae ( ) |
||
real(kind=wp), | intent(in), | dimension(:) | :: | z |
|
|
integer(kind=ip), | intent(in) | :: | nz |
number of abcissae ( ) |
||
real(kind=wp), | intent(in), | dimension(:) | :: | q |
|
|
integer(kind=ip), | intent(in) | :: | nq |
number of abcissae ( ) |
||
real(kind=wp), | intent(in), | dimension(:) | :: | r |
|
|
integer(kind=ip), | intent(in) | :: | nr |
number of abcissae ( ) |
||
real(kind=wp), | intent(in), | dimension(:,:,:,:,:) | :: | fcn |
|
|
integer(kind=ip), | intent(in) | :: | kx |
the order of spline pieces in ( ). (order = polynomial degree + 1) |
||
integer(kind=ip), | intent(in) | :: | ky |
the order of spline pieces in ( ). (order = polynomial degree + 1) |
||
integer(kind=ip), | intent(in) | :: | kz |
the order of spline pieces in ( ). (order = polynomial degree + 1) |
||
integer(kind=ip), | intent(in) | :: | kq |
the order of spline pieces in ( ). (order = polynomial degree + 1) |
||
integer(kind=ip), | intent(in) | :: | kr |
the order of spline pieces in ( ). (order = polynomial degree + 1) |
||
integer(kind=ip), | intent(in) | :: | iknot |
knot sequence flag:
|
||
real(kind=wp), | intent(inout), | dimension(:) | :: | tx |
The
Must be non-decreasing. |
|
real(kind=wp), | intent(inout), | dimension(:) | :: | ty |
The
Must be non-decreasing. |
|
real(kind=wp), | intent(inout), | dimension(:) | :: | tz |
The
Must be non-decreasing. |
|
real(kind=wp), | intent(inout), | dimension(:) | :: | tq |
The
Must be non-decreasing. |
|
real(kind=wp), | intent(inout), | dimension(:) | :: | tr |
The
Must be non-decreasing. |
|
real(kind=wp), | intent(out), | dimension(:,:,:,:,:) | :: | bcoef |
|
|
integer(kind=ip), | intent(out) | :: | iflag |
|
pure subroutine db5ink(x,nx,y,ny,z,nz,q,nq,r,nr,& fcn,& kx,ky,kz,kq,kr,& iknot,& tx,ty,tz,tq,tr,& bcoef,iflag) implicit none integer(ip),intent(in) :: nx !! number of \(x\) abcissae ( \( \ge 3 \) ) integer(ip),intent(in) :: ny !! number of \(y\) abcissae ( \( \ge 3 \) ) integer(ip),intent(in) :: nz !! number of \(z\) abcissae ( \( \ge 3 \) ) integer(ip),intent(in) :: nq !! number of \(q\) abcissae ( \( \ge 3 \) ) integer(ip),intent(in) :: nr !! number of \(r\) abcissae ( \( \ge 3 \) ) integer(ip),intent(in) :: kx !! the order of spline pieces in \(x\) !! ( \( 2 \le k_x < n_x \) ). !! (order = polynomial degree + 1) integer(ip),intent(in) :: ky !! the order of spline pieces in \(y\) !! ( \( 2 \le k_y < n_y \) ). !! (order = polynomial degree + 1) integer(ip),intent(in) :: kz !! the order of spline pieces in \(z\) !! ( \( 2 \le k_z < n_z \) ). !! (order = polynomial degree + 1) integer(ip),intent(in) :: kq !! the order of spline pieces in \(q\) !! ( \( 2 \le k_q < n_q \) ). !! (order = polynomial degree + 1) integer(ip),intent(in) :: kr !! the order of spline pieces in \(r\) !! ( \( 2 \le k_r < n_r \) ). !! (order = polynomial degree + 1) real(wp),dimension(:),intent(in) :: x !! `(nx)` array of \(x\) abcissae. must be strictly increasing. real(wp),dimension(:),intent(in) :: y !! `(ny)` array of \(y\) abcissae. must be strictly increasing. real(wp),dimension(:),intent(in) :: z !! `(nz)` array of \(z\) abcissae. must be strictly increasing. real(wp),dimension(:),intent(in) :: q !! `(nq)` array of \(q\) abcissae. must be strictly increasing. real(wp),dimension(:),intent(in) :: r !! `(nr)` array of \(r\) abcissae. must be strictly increasing. real(wp),dimension(:,:,:,:,:),intent(in) :: fcn !! `(nx,ny,nz,nq,nr)` matrix of function values to interpolate. !! `fcn(i,j,k,q,r)` should contain the function value at the !! point (`x(i)`,`y(j)`,`z(k)`,`q(l)`,`r(m)`) integer(ip),intent(in) :: iknot !! knot sequence flag: !! !! * 0 = knot sequence chosen by [[db5ink]]. !! * 1 = knot sequence chosen by user. real(wp),dimension(:),intent(inout) :: tx !! The `(nx+kx)` knots in the \(x\) direction for the spline !! interpolant. !! !! * If `iknot=0` these are chosen by [[db5ink]]. !! * If `iknot=1` these are specified by the user. !! !! Must be non-decreasing. real(wp),dimension(:),intent(inout) :: ty !! The `(ny+ky)` knots in the \(y\) direction for the spline !! interpolant. !! !! * If `iknot=0` these are chosen by [[db5ink]]. !! * If `iknot=1` these are specified by the user. !! !! Must be non-decreasing. real(wp),dimension(:),intent(inout) :: tz !! The `(nz+kz)` knots in the \(z\) direction for the spline !! interpolant. !! !! * If `iknot=0` these are chosen by [[db5ink]]. !! * If `iknot=1` these are specified by the user. !! !! Must be non-decreasing. real(wp),dimension(:),intent(inout) :: tq !! The `(nq+kq)` knots in the \(q\) direction for the spline !! interpolant. !! !! * If `iknot=0` these are chosen by [[db5ink]]. !! * If `iknot=1` these are specified by the user. !! !! Must be non-decreasing. real(wp),dimension(:),intent(inout) :: tr !! The `(nr+kr)` knots in the \(r\) direction for the spline !! interpolant. !! !! * If `iknot=0` these are chosen by [[db5ink]]. !! * If `iknot=1` these are specified by the user. !! !! Must be non-decreasing. real(wp),dimension(:,:,:,:,:),intent(out) :: bcoef !! `(nx,ny,nz,nq,nr)` matrix of coefficients of the b-spline !! interpolant. integer(ip),intent(out) :: iflag !! * 0 = successful execution. !! * 2 = `iknot` out of range. !! * 3 = `nx` out of range. !! * 4 = `kx` out of range. !! * 5 = `x` not strictly increasing. !! * 6 = `tx` not non-decreasing. !! * 7 = `ny` out of range. !! * 8 = `ky` out of range. !! * 9 = `y` not strictly increasing. !! * 10 = `ty` not non-decreasing. !! * 11 = `nz` out of range. !! * 12 = `kz` out of range. !! * 13 = `z` not strictly increasing. !! * 14 = `tz` not non-decreasing. !! * 15 = `nq` out of range. !! * 16 = `kq` out of range. !! * 17 = `q` not strictly increasing. !! * 18 = `tq` not non-decreasing. !! * 19 = `nr` out of range. !! * 20 = `kr` out of range. !! * 21 = `r` not strictly increasing. !! * 22 = `tr` not non-decreasing. !! * 700 = `size(x)` \( \ne \) `size(fcn,1)` !! * 701 = `size(y)` \( \ne \) `size(fcn,2)` !! * 702 = `size(z)` \( \ne \) `size(fcn,3)` !! * 703 = `size(q)` \( \ne \) `size(fcn,4)` !! * 704 = `size(r)` \( \ne \) `size(fcn,5)` !! * 706 = `size(x)` \( \ne \) `nx` !! * 707 = `size(y)` \( \ne \) `ny` !! * 708 = `size(z)` \( \ne \) `nz` !! * 709 = `size(q)` \( \ne \) `nq` !! * 710 = `size(r)` \( \ne \) `nr` !! * 712 = `size(tx)` \( \ne \) `nx+kx` !! * 713 = `size(ty)` \( \ne \) `ny+ky` !! * 714 = `size(tz)` \( \ne \) `nz+kz` !! * 715 = `size(tq)` \( \ne \) `nq+kq` !! * 716 = `size(tr)` \( \ne \) `nr+kr` !! * 800 = `size(x)` \( \ne \) `size(bcoef,1)` !! * 801 = `size(y)` \( \ne \) `size(bcoef,2)` !! * 802 = `size(z)` \( \ne \) `size(bcoef,3)` !! * 803 = `size(q)` \( \ne \) `size(bcoef,4)` !! * 804 = `size(r)` \( \ne \) `size(bcoef,5)` logical :: status_ok real(wp),dimension(:),allocatable :: temp !! work array of length `nx*ny*nz*nq*nr` real(wp),dimension(:),allocatable :: work !! work array of length `max(2*kx*(nx+1), !! 2*ky*(ny+1),2*kz*(nz+1),2*kq*(nq+1),2*kr*(nr+1))` integer(ip) :: i, j, k, l, m, ii !! counter ! check validity of input call check_inputs( iknot,& iflag,& nx=nx,ny=ny,nz=nz,nq=nq,nr=nr,& kx=kx,ky=ky,kz=kz,kq=kq,kr=kr,& x=x,y=y,z=z,q=q,r=r,& tx=tx,ty=ty,tz=tz,tq=tq,tr=tr,& f5=fcn,& bcoef5=bcoef,& status_ok=status_ok) if (status_ok) then ! choose knots if (iknot == 0_ip) then call dbknot(x,nx,kx,tx) call dbknot(y,ny,ky,ty) call dbknot(z,nz,kz,tz) call dbknot(q,nq,kq,tq) call dbknot(r,nr,kr,tr) end if allocate(temp(nx*ny*nz*nq*nr)) allocate(work(max(2_ip*kx*(nx+1_ip),2_ip*ky*(ny+1_ip),2_ip*kz*(nz+1_ip),2_ip*kq*(nq+1_ip),2_ip*kr*(nr+1_ip)))) ! copy fcn to work in packed for dbtpcf !temp(1:nx*ny*nz*nq*nr) = reshape( fcn, [nx*ny*nz*nq*nr] ) ! replaced with loops to avoid stack ! overflow for large data set: ii = 0_ip do m = 1_ip, nr do l = 1_ip, nq do k = 1_ip, nz do j = 1_ip, ny do i = 1_ip, nx ii = ii + 1_ip temp(ii) = fcn(i,j,k,l,m) end do end do end do end do end do ! construct b-spline coefficients call dbtpcf(x,nx,temp, nx,ny*nz*nq*nr,tx,kx,bcoef,work,iflag) if (iflag==0_ip) call dbtpcf(y,ny,bcoef, ny,nx*nz*nq*nr,ty,ky,temp, work,iflag) if (iflag==0_ip) call dbtpcf(z,nz,temp, nz,nx*ny*nq*nr,tz,kz,bcoef,work,iflag) if (iflag==0_ip) call dbtpcf(q,nq,bcoef, nq,nx*ny*nz*nr,tq,kq,temp, work,iflag) if (iflag==0_ip) call dbtpcf(r,nr,temp, nr,nx*ny*nz*nq,tr,kr,bcoef,work,iflag) deallocate(temp) deallocate(work) end if end subroutine db5ink