46 |
|
!! PURPOSE: Generic Spline interpolation routines. |
47 |
|
!! |
48 |
|
!! @author Charles F. Vardeman II |
49 |
< |
!! @version $Id: interpolation.F90,v 1.7 2006-04-20 18:24:24 gezelter Exp $ |
49 |
> |
!! @version $Id: interpolation.F90,v 1.9 2006-06-06 17:43:28 gezelter Exp $ |
50 |
|
|
51 |
|
|
52 |
|
module interpolation |
224 |
|
cs%b(n) = cs%b(n-1) + h(n-1) * (2.0_dp*cs%c(n-1) + h(n-1)*3.0_dp*cs%d(n-1)) |
225 |
|
|
226 |
|
if (isUniform) then |
227 |
< |
cs%dx_i = 1.0d0 / (x(2) - x(1)) |
227 |
> |
cs%dx_i = 1.0_dp / (x(2) - x(1)) |
228 |
|
endif |
229 |
|
|
230 |
|
return |
252 |
|
|
253 |
|
type(cubicSpline) :: this |
254 |
|
|
255 |
< |
if(associated(this%x)) then |
256 |
< |
deallocate(this%x) |
257 |
< |
this%x => null() |
255 |
> |
if(associated(this%d)) then |
256 |
> |
deallocate(this%d) |
257 |
> |
this%d => null() |
258 |
|
end if |
259 |
|
if(associated(this%c)) then |
260 |
|
deallocate(this%c) |
261 |
|
this%c => null() |
262 |
+ |
end if |
263 |
+ |
if(associated(this%b)) then |
264 |
+ |
deallocate(this%b) |
265 |
+ |
this%b => null() |
266 |
+ |
end if |
267 |
+ |
if(associated(this%y)) then |
268 |
+ |
deallocate(this%y) |
269 |
+ |
this%y => null() |
270 |
|
end if |
271 |
+ |
if(associated(this%x)) then |
272 |
+ |
deallocate(this%x) |
273 |
+ |
this%x => null() |
274 |
+ |
end if |
275 |
|
|
276 |
|
this%n = 0 |
277 |
|
|
322 |
|
! Find the interval J = [ cs%x(J), cs%x(J+1) ] that contains |
323 |
|
! or is nearest to xval. |
324 |
|
|
325 |
< |
j = MAX(1, MIN(cs%n-1, idint((xval-cs%x(1)) * cs%dx_i) + 1)) |
325 |
> |
j = MAX(1, MIN(cs%n-1, int((xval-cs%x(1)) * cs%dx_i) + 1)) |
326 |
|
|
327 |
|
dx = xval - cs%x(j) |
328 |
|
yval = cs%y(j) + dx*(cs%b(j) + dx*(cs%c(j) + dx*cs%d(j))) |
344 |
|
! or is nearest to xval. |
345 |
|
|
346 |
|
|
347 |
< |
j = MAX(1, MIN(cs%n-1, idint((xval-cs%x(1)) * cs%dx_i) + 1)) |
347 |
> |
j = MAX(1, MIN(cs%n-1, int((xval-cs%x(1)) * cs%dx_i) + 1)) |
348 |
|
|
349 |
|
dx = xval - cs%x(j) |
350 |
|
yval = cs%y(j) + dx*(cs%b(j) + dx*(cs%c(j) + dx*cs%d(j))) |
351 |
|
|
352 |
< |
dydx = cs%b(j) + dx*(2.0d0 * cs%c(j) + 3.0d0 * dx * cs%d(j)) |
352 |
> |
dydx = cs%b(j) + dx*(2.0_dp * cs%c(j) + 3.0_dp * dx * cs%d(j)) |
353 |
|
|
354 |
|
return |
355 |
|
end subroutine lookupUniformSpline1d |