V3FIT
ezspline_gradient.f90
1 !/////
2 ! R8 !
3 !/////
4 !
5 ! 1-D
6 !
7 
8 subroutine ezspline_gradient1_r8(spline_o, &
9  p1, df, ier)
10  use ezspline_obj
11  implicit none
12  type(EZspline1_r8) spline_o
13  real(ezspline_r8), intent(in) :: p1
14  real(ezspline_r8), intent(out) :: df
15  integer, intent(out) :: ier
16 
17  integer ifail
18  integer, parameter :: ict(3) = (/0, 1, 0/)
19 
20  ier = 0
21  ifail=0
22  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
23  ier = 94
24  return
25  endif
26 
27  if(spline_o%isLinear == 1) then
28 
29  call r8pc1ev(p1, &
30  & spline_o%x1(1), spline_o%n1, &
31  & spline_o%ilin1, &
32  & spline_o%fspl(1,1), &
33  & ict, df, ifail)
34 
35  else if(spline_o%isHermite == 0) then
36 
37  call r8evspline(p1, &
38  & spline_o%x1(1), spline_o%n1, &
39  & spline_o%ilin1, &
40  & spline_o%fspl(1,1), &
41  & ict, df, ifail)
42 
43  else
44 
45  call r8herm1ev(p1, &
46  & spline_o%x1(1), spline_o%n1, &
47  & spline_o%ilin1, &
48  & spline_o%fspl(1,1), &
49  & ict, df, ifail)
50 
51  endif
52  if(ifail/=0) ier = 95
53 
54 end subroutine ezspline_gradient1_r8
55 
56 
57 subroutine ezspline_gradient1_array_r8(spline_o, k1, &
58  p1, df, ier)
59  use ezspline_obj
60  implicit none
61  type(EZspline1_r8) spline_o
62  integer, intent(in) :: k1
63  real(ezspline_r8), intent(in) :: p1(k1)
64  real(ezspline_r8), intent(out) :: df(k1)
65  !
66  ! ier=95 some error occurred in EZspline_gradient
67  !
68  integer, intent(out) :: ier
69  integer ifail
70  integer, parameter :: ict(3) = (/0, 1, 0/)
71  integer:: iwarn=0
72 
73  ier = 0
74  ifail=0
75  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
76  ier = 94
77  return
78  endif
79 
80  if(spline_o%isLinear == 1) then
81 
82  call r8vecpc1(ict, k1, p1, k1, df, &
83  & spline_o%n1, spline_o%x1pkg(1,1), &
84  & spline_o%fspl(1,1), &
85  & iwarn, ifail)
86 
87  else if(spline_o%isHermite == 0) then
88 
89  call r8vecspline(ict, k1, p1, k1, df, &
90  & spline_o%n1,spline_o%x1pkg(1,1), &
91  & spline_o%fspl(1,1), &
92  & iwarn, ifail)
93 
94 
95  else
96 
97  call r8vecherm1(ict, k1, p1, k1, df, &
98  & spline_o%n1, spline_o%x1pkg(1,1), &
99  & spline_o%fspl(1,1), &
100  & iwarn, ifail)
101 
102  endif
103 
104  if(ifail /= 0) ier = 95
105 end subroutine ezspline_gradient1_array_r8
106 
107 
108 !!
109 !! 2-D
110 !!
111 
112 subroutine ezspline_gradient2_r8(spline_o, &
113  p1, p2, df, ier)
114  use ezspline_obj
115  implicit none
116  type(EZspline2_r8) spline_o
117  real(ezspline_r8), intent(in) :: p1, p2
118  real(ezspline_r8), intent(out) :: df(2)
119  !
120  integer, intent(out) :: ier
121  integer ifail
122  integer, parameter :: ict(6) = (/0, 1, 1, 0, 0, 0/)
123 
124  ier = 0
125  ifail=0
126  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
127  ier = 94
128  return
129  endif
130 
131  if (spline_o%isHybrid == 1) then
132 
133  call r8evintrp2d(p1, p2, &
134  & spline_o%x1(1), spline_o%n1, &
135  & spline_o%x2(1), spline_o%n2, &
136  & spline_o%hspline, spline_o%fspl(1,1,1), &
137  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
138  & size(spline_o%fspl,3), &
139  & ict, df, ifail)
140 
141  else if(spline_o%isLinear == 1) then
142 
143  call r8pc2ev(p1, p2, &
144  & spline_o%x1(1), spline_o%n1, &
145  & spline_o%x2(1), spline_o%n2, &
146  & spline_o%ilin1, spline_o%ilin2, &
147  & spline_o%fspl(1,1,1), spline_o%n1, &
148  & ict, df, ifail)
149 
150  else if(spline_o%isHermite == 0) then
151 
152  call r8evbicub(p1, p2, &
153  & spline_o%x1(1), spline_o%n1, &
154  & spline_o%x2(1), spline_o%n2, &
155  & spline_o%ilin1, spline_o%ilin2, &
156  & spline_o%fspl(1,1,1), spline_o%n1, &
157  & ict, df, ifail)
158 
159  else
160 
161  call r8herm2ev(p1, p2, &
162  & spline_o%x1(1), spline_o%n1, &
163  & spline_o%x2(1), spline_o%n2, &
164  & spline_o%ilin1, spline_o%ilin2, &
165  & spline_o%fspl(1,1,1), spline_o%n1, &
166  & ict, df, ifail)
167  endif
168 
169  if(ifail/=0) ier = 95
170 
171 end subroutine ezspline_gradient2_r8
172 
173 subroutine ezspline_gradient2_cloud_r8(spline_o, k, &
174  p1, p2, df, ier)
175  use ezspline_obj
176  implicit none
177  type(EZspline2_r8) spline_o
178  integer, intent(in) :: k
179  real(ezspline_r8), intent(in) :: p1(k), p2(k)
180  real(ezspline_r8), intent(out) :: df(k,2)
181  !
182  integer, intent(out) :: ier
183  integer ifail
184  integer, parameter :: ict(6) = (/0, 1, 1, 0, 0, 0/)
185  integer:: iwarn=0
186 
187  ier = 0
188  ifail=0
189  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
190  ier = 94
191  return
192  endif
193 
194  if (spline_o%isHybrid == 1) then
195 
196  call r8vecintrp2d(ict, k, p1, p2, k, df, &
197  & spline_o%n1, spline_o%x1pkg(1,1), &
198  & spline_o%n2, spline_o%x2pkg(1,1), &
199  & spline_o%hspline, spline_o%fspl(1,1,1), &
200  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
201  & size(spline_o%fspl,3), &
202  & iwarn, ifail)
203 
204  else if(spline_o%isLinear == 1) then
205 
206  call r8vecpc2(ict, k, p1, p2, k, df, &
207  & spline_o%n1, spline_o%x1pkg(1,1), &
208  & spline_o%n2, spline_o%x2pkg(1,1), &
209  & spline_o%fspl(1,1,1), spline_o%n1, &
210  & iwarn,ifail)
211 
212  else if(spline_o%isHermite == 0) then
213 
214  call r8vecbicub(ict, k, p1, p2, k, df, &
215  & spline_o%n1,spline_o%x1pkg(1,1), &
216  & spline_o%n2,spline_o%x2pkg(1,1), &
217  & spline_o%fspl(1,1,1), spline_o%n1, &
218  & iwarn, ifail)
219 
220  else
221 
222  call r8vecherm2(ict, k, p1, p2, k, df, &
223  & spline_o%n1, spline_o%x1pkg(1,1), &
224  & spline_o%n2, spline_o%x2pkg(1,1), &
225  & spline_o%fspl(1,1,1), spline_o%n1, &
226  & iwarn,ifail)
227 
228  endif
229  if(ifail /= 0) ier = 95
230 
231 end subroutine ezspline_gradient2_cloud_r8
232 
233 
234 subroutine ezspline_gradient2_array_r8(spline_o, k1, k2, &
235  p1, p2, df, ier)
236  use ezspline_obj
237  implicit none
238  type(EZspline2_r8) spline_o
239  integer, intent(in) :: k1, k2
240  real(ezspline_r8), intent(in) :: p1(k1), p2(k2)
241  real(ezspline_r8), intent(out) :: df(k1, k2, 2)
242  !
243  integer, intent(out) :: ier
244  integer ifail
245  integer, parameter :: ict(6) = (/0, 1, 1, 0, 0, 0/)
246  integer:: iwarn=0
247  real(ezspline_r8), dimension(:), allocatable :: p1_cloud, p2_cloud
248  integer k12
249 
250  ier = 0
251  ifail=0
252  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
253  ier = 94
254  return
255  endif
256 
257 
258  k12 = k1*k2
259  allocate(p1_cloud(k12), p2_cloud(k12), stat=ifail)
260  if(ifail /= 0) then
261  ier = 30
262  return
263  endif
264 
265  p1_cloud = reshape( &
266  & source=spread(source=p1, dim=2, ncopies=k2), shape=(/k12/))
267  p2_cloud = reshape( &
268  & source=spread(source=p2, dim=1, ncopies=k1), shape=(/k12/))
269 
270  if(spline_o%isHybrid == 1) then
271 
272  call r8vecintrp2d(ict, k12, p1_cloud, p2_cloud, k12, df, &
273  & spline_o%n1, spline_o%x1pkg(1,1), &
274  & spline_o%n2, spline_o%x2pkg(1,1), &
275  & spline_o%hspline, spline_o%fspl(1,1,1), &
276  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
277  & size(spline_o%fspl,3), &
278  & iwarn, ifail)
279 
280  else if(spline_o%isLinear == 1) then
281 
282  call r8vecpc2(ict, k12, p1_cloud, p2_cloud, k12, df, &
283  & spline_o%n1, spline_o%x1pkg(1,1), &
284  & spline_o%n2, spline_o%x2pkg(1,1), &
285  & spline_o%fspl(1,1,1), spline_o%n1, &
286  & iwarn, ifail)
287 
288  else if(spline_o%isHermite == 0) then
289 
290  call r8vecbicub(ict, k12, p1_cloud, p2_cloud, k12, df, &
291  & spline_o%n1,spline_o%x1pkg(1,1), &
292  & spline_o%n2,spline_o%x2pkg(1,1), &
293  & spline_o%fspl(1,1,1), spline_o%n1, &
294  & iwarn, ifail)
295 
296  else
297 
298  call r8vecherm2(ict, k12, p1_cloud, p2_cloud, k12, df, &
299  & spline_o%n1, spline_o%x1pkg(1,1), &
300  & spline_o%n2, spline_o%x2pkg(1,1), &
301  & spline_o%fspl(1,1,1), spline_o%n1, &
302  & iwarn, ifail)
303 
304  endif
305 
306  deallocate(p1_cloud, p2_cloud, stat=ifail)
307  if(ifail /= 0) then
308  ier = 31
309  return
310  endif
311 
312  if(ifail /= 0) ier = 95
313 
314 end subroutine ezspline_gradient2_array_r8
315 
316 
317 !!!
318 !!! 3-D
319 !!!
320 
321 subroutine ezspline_gradient3_r8(spline_o, &
322  p1, p2, p3, df, ier)
323  use ezspline_obj
324  implicit none
325  type(EZspline3_r8) spline_o
326  real(ezspline_r8), intent(in) :: p1, p2, p3
327  real(ezspline_r8), intent(out) :: df(3)
328  !
329  integer, intent(out) :: ier
330  integer ifail
331  integer, parameter :: ict(10) = (/0, 1, 1, 1, 0, 0, 0, 0, 0, 0/)
332 
333  ier = 0
334  ifail=0
335  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
336  ier = 94
337  return
338  endif
339 
340  if (spline_o%isHybrid == 1) then
341 
342  call r8evintrp3d(p1, p2, p3, &
343  & spline_o%x1(1), spline_o%n1, &
344  & spline_o%x2(1), spline_o%n2, &
345  & spline_o%x3(1), spline_o%n3, &
346  & spline_o%hspline, spline_o%fspl(1,1,1,1), &
347  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
348  & size(spline_o%fspl,3), size(spline_o%fspl,4), &
349  & ict, df, ifail)
350 
351  else if(spline_o%isLinear == 1) then
352 
353  call r8pc3ev(p1, p2, p3, &
354  & spline_o%x1(1), spline_o%n1, &
355  & spline_o%x2(1), spline_o%n2, &
356  & spline_o%x3(1), spline_o%n3, &
357  & spline_o%ilin1, spline_o%ilin2, spline_o%ilin3, &
358  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
359  & ict, df, ifail)
360 
361  else if(spline_o%isHermite == 0) then
362 
363  call r8evtricub(p1, p2, p3, &
364  & spline_o%x1(1), spline_o%n1, &
365  & spline_o%x2(1), spline_o%n2, &
366  & spline_o%x3(1), spline_o%n3, &
367  & spline_o%ilin1, spline_o%ilin2, spline_o%ilin3, &
368  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
369  & ict, df, ifail)
370 
371  else
372 
373  call r8herm3ev(p1, p2, p3, &
374  & spline_o%x1(1), spline_o%n1, &
375  & spline_o%x2(1), spline_o%n2, &
376  & spline_o%x3(1), spline_o%n3, &
377  & spline_o%ilin1, spline_o%ilin2, spline_o%ilin3, &
378  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
379  & ict, df, ifail)
380 
381  endif
382 
383  if(ifail/=0) ier = 95
384 
385 end subroutine ezspline_gradient3_r8
386 
387 subroutine ezspline_gradient3_cloud_r8(spline_o, k, &
388  p1, p2, p3, df, ier)
389  use ezspline_obj
390  implicit none
391  type(EZspline3_r8) spline_o
392  integer, intent(in) :: k
393  real(ezspline_r8), intent(in) :: p1(k), p2(k), p3(k)
394  real(ezspline_r8), intent(out) :: df(k, 3)
395  !
396  integer, intent(out) :: ier
397  integer ifail
398  integer, parameter :: ict(10) = (/0, 1, 1, 1, 0, 0, 0, 0, 0, 0/)
399  integer :: iwarn=0
400 
401  ier = 0
402  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
403  ier = 94
404  return
405  endif
406 
407  if (spline_o%isHybrid == 1) then
408 
409  call r8vecintrp3d(ict, k, p1, p2, p3, k, df, &
410  & spline_o%n1, spline_o%x1pkg(1,1), &
411  & spline_o%n2, spline_o%x2pkg(1,1), &
412  & spline_o%n3, spline_o%x3pkg(1,1), &
413  & spline_o%hspline, spline_o%fspl(1,1,1,1), &
414  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
415  & size(spline_o%fspl,3), size(spline_o%fspl,4), &
416  & iwarn, ifail)
417 
418  else if(spline_o%isLinear == 1) then
419 
420  call r8vecpc3(ict, k, p1, p2, p3, k, df, &
421  & spline_o%n1, spline_o%x1pkg(1,1), &
422  & spline_o%n2, spline_o%x2pkg(1,1), &
423  & spline_o%n3, spline_o%x3pkg(1,1), &
424  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
425  & iwarn,ifail)
426 
427  else if(spline_o%isHermite == 0) then
428 
429  call r8vectricub(ict, k, p1, p2, p3, k, df, &
430  & spline_o%n1,spline_o%x1pkg(1,1), &
431  & spline_o%n2,spline_o%x2pkg(1,1), &
432  & spline_o%n3,spline_o%x3pkg(1,1), &
433  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
434  & iwarn, ifail)
435 
436  else
437 
438  call r8vecherm3(ict, k, p1, p2, p3, k, df, &
439  & spline_o%n1, spline_o%x1pkg(1,1), &
440  & spline_o%n2, spline_o%x2pkg(1,1), &
441  & spline_o%n3, spline_o%x3pkg(1,1), &
442  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
443  & iwarn,ifail)
444 
445  endif
446 
447  if(ifail /= 0) ier = 95
448 
449 
450 end subroutine ezspline_gradient3_cloud_r8
451 
452 
453 subroutine ezspline_gradient3_array_r8(spline_o, k1, k2, k3, &
454  p1, p2, p3, df, ier)
455  use ezspline_obj
456  implicit none
457  type(EZspline3_r8) spline_o
458  integer, intent(in) :: k1, k2, k3
459  real(ezspline_r8), intent(in) :: p1(k1), p2(k2), p3(k3)
460  real(ezspline_r8), intent(out) :: df(k1, k2, k3, 3)
461  !
462  integer, intent(out) :: ier
463  integer ifail
464  integer, parameter :: ict(10) = (/0, 1, 1, 1, 0, 0, 0, 0, 0, 0/)
465  integer:: iwarn=0
466  real(ezspline_r8), dimension(:), allocatable :: p1_cloud, p2_cloud, p3_cloud
467  integer k123
468 
469  ier = 0
470  ifail=0
471  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
472  ier = 94
473  return
474  endif
475 
476  k123 = k1*k2*k3
477  allocate(p1_cloud(k123), p2_cloud(k123), p3_cloud(k123), stat=ifail)
478  if(ifail /= 0) then
479  ier = 30
480  return
481  endif
482 
483  p1_cloud = reshape(source=spread( &
484  & source=spread(source=p1, dim=2, ncopies=k2), &
485  & dim=3, ncopies=k3), shape=(/k123/))
486  p2_cloud = reshape(source=spread( &
487  & source=spread(source=p2, dim=1, ncopies=k1), &
488  & dim=3, ncopies=k3), shape=(/k123/))
489  p3_cloud = reshape(source=spread( &
490  & source=spread(source=p3, dim=1, ncopies=k1), &
491  & dim=2, ncopies=k2), shape=(/k123/))
492 
493  if (spline_o%isHybrid == 1) then
494 
495  call r8vecintrp3d(ict, k123, p1_cloud, p2_cloud, p3_cloud, k123, df, &
496  & spline_o%n1, spline_o%x1pkg(1,1), &
497  & spline_o%n2, spline_o%x2pkg(1,1), &
498  & spline_o%n3, spline_o%x3pkg(1,1), &
499  & spline_o%hspline, spline_o%fspl(1,1,1,1), &
500  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
501  & size(spline_o%fspl,3), size(spline_o%fspl,4), &
502  & iwarn, ifail)
503 
504  else if(spline_o%isLinear == 1) then
505 
506  call r8vecpc3(ict, k123, p1_cloud, p2_cloud, p3_cloud, k123, df, &
507  & spline_o%n1, spline_o%x1pkg(1,1), &
508  & spline_o%n2, spline_o%x2pkg(1,1), &
509  & spline_o%n3, spline_o%x3pkg(1,1), &
510  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
511  & iwarn,ifail)
512 
513  else if(spline_o%isHermite == 0) then
514 
515  call r8vectricub(ict, k123, p1_cloud, p2_cloud, p3_cloud, k123, df, &
516  & spline_o%n1,spline_o%x1pkg(1,1), &
517  & spline_o%n2,spline_o%x2pkg(1,1), &
518  & spline_o%n3,spline_o%x3pkg(1,1), &
519  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
520  & iwarn, ifail)
521 
522  else
523 
524  call r8vecherm3(ict, k123, p1_cloud, p2_cloud, p3_cloud, k123, df, &
525  & spline_o%n1, spline_o%x1pkg(1,1), &
526  & spline_o%n2, spline_o%x2pkg(1,1), &
527  & spline_o%n3, spline_o%x3pkg(1,1), &
528  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
529  & iwarn,ifail)
530 
531  endif
532 
533  if(ifail /= 0) ier = 95
534 
535  deallocate(p1_cloud, p2_cloud, p3_cloud, stat=ifail)
536  if(ifail /= 0) then
537  ier = 31
538  return
539  endif
540 
541 end subroutine ezspline_gradient3_array_r8
542 !/////
543 ! R4 !
544 !/////
545 !
546 ! 1-D
547 !
548 
549 
550 subroutine ezspline_gradient1_r4(spline_o, &
551  p1, df, ier)
552  use ezspline_obj
553  implicit none
554  type(EZspline1_r4) spline_o
555  real(ezspline_r4), intent(in) :: p1
556  real(ezspline_r4), intent(out) :: df
557  integer, intent(out) :: ier
558 
559  integer ifail
560  integer, parameter :: ict(3) = (/0, 1, 0/)
561 
562  ier = 0
563  ifail=0
564  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
565  ier = 94
566  return
567  endif
568 
569  if(spline_o%isLinear == 1) then
570 
571  call pc1ev(p1, &
572  & spline_o%x1(1), spline_o%n1, &
573  & spline_o%ilin1, &
574  & spline_o%fspl(1,1), &
575  & ict, df, ifail)
576 
577  else if(spline_o%isHermite == 0) then
578 
579  call evspline(p1, &
580  & spline_o%x1(1), spline_o%n1, &
581  & spline_o%ilin1, &
582  & spline_o%fspl(1,1), &
583  & ict, df, ifail)
584 
585  else
586 
587  call herm1ev(p1, &
588  & spline_o%x1(1), spline_o%n1, &
589  & spline_o%ilin1, &
590  & spline_o%fspl(1,1), &
591  & ict, df, ifail)
592 
593  endif
594  if(ifail/=0) ier = 95
595 
596 end subroutine ezspline_gradient1_r4
597 
598 
599 subroutine ezspline_gradient1_array_r4(spline_o, k1, &
600  p1, df, ier)
601  use ezspline_obj
602  implicit none
603  type(EZspline1_r4) spline_o
604  integer, intent(in) :: k1
605  real(ezspline_r4), intent(in) :: p1(k1)
606  real(ezspline_r4), intent(out) :: df(k1)
607  !
608  ! ier=95 some error occurred in EZspline_gradient
609  !
610  integer, intent(out) :: ier
611  integer ifail
612  integer, parameter :: ict(3) = (/0, 1, 0/)
613  integer:: iwarn=0
614 
615  ier = 0
616  ifail=0
617  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
618  ier = 94
619  return
620  endif
621 
622  if(spline_o%isLinear == 1) then
623 
624  call vecpc1(ict, k1, p1, k1, df, &
625  & spline_o%n1, spline_o%x1pkg(1,1), &
626  & spline_o%fspl(1,1), &
627  & iwarn, ifail)
628 
629  else if(spline_o%isHermite == 0) then
630 
631  call vecspline(ict, k1, p1, k1, df, &
632  & spline_o%n1,spline_o%x1pkg(1,1), &
633  & spline_o%fspl(1,1), &
634  & iwarn, ifail)
635 
636 
637  else
638 
639  call vecherm1(ict, k1, p1, k1, df, &
640  & spline_o%n1, spline_o%x1pkg(1,1), &
641  & spline_o%fspl(1,1), &
642  & iwarn, ifail)
643 
644  endif
645 
646  if(ifail /= 0) ier = 95
647 end subroutine ezspline_gradient1_array_r4
648 
649 
650 !!
651 !! 2-D
652 !!
653 
654 subroutine ezspline_gradient2_r4(spline_o, &
655  p1, p2, df, ier)
656  use ezspline_obj
657  implicit none
658  type(EZspline2_r4) spline_o
659  real(ezspline_r4), intent(in) :: p1, p2
660  real(ezspline_r4), intent(out) :: df(2)
661  !
662  integer, intent(out) :: ier
663  integer ifail
664  integer, parameter :: ict(6) = (/0, 1, 1, 0, 0, 0/)
665 
666  ier = 0
667  ifail=0
668  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
669  ier = 94
670  return
671  endif
672 
673  if (spline_o%isHybrid == 1) then
674 
675  call evintrp2d(p1, p2, &
676  & spline_o%x1(1), spline_o%n1, &
677  & spline_o%x2(1), spline_o%n2, &
678  & spline_o%hspline, spline_o%fspl(1,1,1), &
679  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
680  & size(spline_o%fspl,3), &
681  & ict, df, ifail)
682 
683  else if(spline_o%isLinear == 1) then
684 
685  call pc2ev(p1, p2, &
686  & spline_o%x1(1), spline_o%n1, &
687  & spline_o%x2(1), spline_o%n2, &
688  & spline_o%ilin1, spline_o%ilin2, &
689  & spline_o%fspl(1,1,1), spline_o%n1, &
690  & ict, df, ifail)
691 
692  else if(spline_o%isHermite == 0) then
693 
694  call evbicub(p1, p2, &
695  & spline_o%x1(1), spline_o%n1, &
696  & spline_o%x2(1), spline_o%n2, &
697  & spline_o%ilin1, spline_o%ilin2, &
698  & spline_o%fspl(1,1,1), spline_o%n1, &
699  & ict, df, ifail)
700 
701  else
702 
703  call herm2ev(p1, p2, &
704  & spline_o%x1(1), spline_o%n1, &
705  & spline_o%x2(1), spline_o%n2, &
706  & spline_o%ilin1, spline_o%ilin2, &
707  & spline_o%fspl(1,1,1), spline_o%n1, &
708  & ict, df, ifail)
709 
710  endif
711 
712  if(ifail/=0) ier = 95
713 
714 end subroutine ezspline_gradient2_r4
715 
716 subroutine ezspline_gradient2_cloud_r4(spline_o, k, &
717  p1, p2, df, ier)
718  use ezspline_obj
719  implicit none
720  type(EZspline2_r4) spline_o
721  integer, intent(in) :: k
722  real(ezspline_r4), intent(in) :: p1(k), p2(k)
723  real(ezspline_r4), intent(out) :: df(k,2)
724  !
725  integer, intent(out) :: ier
726  integer ifail
727  integer, parameter :: ict(6) = (/0, 1, 1, 0, 0, 0/)
728  integer:: iwarn=0
729 
730  ier = 0
731  ifail=0
732  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
733  ier = 94
734  return
735  endif
736 
737 
738  if (spline_o%isHybrid == 1) then
739 
740  call vecintrp2d(ict, k, p1, p2, k, df, &
741  & spline_o%n1, spline_o%x1pkg(1,1), &
742  & spline_o%n2, spline_o%x2pkg(1,1), &
743  & spline_o%hspline, spline_o%fspl(1,1,1), &
744  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
745  & size(spline_o%fspl,3), &
746  & iwarn, ifail)
747 
748  else if(spline_o%isLinear == 1) then
749 
750  call vecpc2(ict, k, p1, p2, k, df, &
751  & spline_o%n1, spline_o%x1pkg(1,1), &
752  & spline_o%n2, spline_o%x2pkg(1,1), &
753  & spline_o%fspl(1,1,1), spline_o%n1, &
754  & iwarn,ifail)
755 
756  else if(spline_o%isHermite == 0) then
757 
758  call vecbicub(ict, k, p1, p2, k, df, &
759  & spline_o%n1,spline_o%x1pkg(1,1), &
760  & spline_o%n2,spline_o%x2pkg(1,1), &
761  & spline_o%fspl(1,1,1), spline_o%n1, &
762  & iwarn, ifail)
763 
764 
765  else
766 
767  call vecherm2(ict, k, p1, p2, k, df, &
768  & spline_o%n1, spline_o%x1pkg(1,1), &
769  & spline_o%n2, spline_o%x2pkg(1,1), &
770  & spline_o%fspl(1,1,1), spline_o%n1, &
771  & iwarn,ifail)
772 
773  endif
774  if(ifail /= 0) ier = 95
775 
776 end subroutine ezspline_gradient2_cloud_r4
777 
778 
779 subroutine ezspline_gradient2_array_r4(spline_o, k1, k2, &
780  p1, p2, df, ier)
781  use ezspline_obj
782  implicit none
783  type(EZspline2_r4) spline_o
784  integer, intent(in) :: k1, k2
785  real(ezspline_r4), intent(in) :: p1(k1), p2(k2)
786  real(ezspline_r4), intent(out) :: df(k1, k2, 2)
787  !
788  integer, intent(out) :: ier
789  integer ifail
790  integer, parameter :: ict(6) = (/0, 1, 1, 0, 0, 0/)
791  integer:: iwarn=0
792  real(ezspline_r4), dimension(:), allocatable :: p1_cloud, p2_cloud
793  integer k12
794 
795  ier = 0
796  ifail=0
797  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
798  ier = 94
799  return
800  endif
801 
802 
803  k12 = k1*k2
804  allocate(p1_cloud(k12), p2_cloud(k12), stat=ifail)
805  if(ifail /= 0) then
806  ier = 30
807  return
808  endif
809 
810  p1_cloud = reshape( &
811  & source=spread(source=p1, dim=2, ncopies=k2), shape=(/k12/))
812  p2_cloud = reshape( &
813  & source=spread(source=p2, dim=1, ncopies=k1), shape=(/k12/))
814 
815  if (spline_o%isHybrid == 1) then
816 
817  call vecintrp2d(ict, k12, p1_cloud, p2_cloud, k12, df, &
818  & spline_o%n1, spline_o%x1pkg(1,1), &
819  & spline_o%n2, spline_o%x2pkg(1,1), &
820  & spline_o%hspline, spline_o%fspl(1,1,1), &
821  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
822  & size(spline_o%fspl,3), &
823  & iwarn, ifail)
824 
825  else if(spline_o%isLinear == 1) then
826 
827  call vecpc2(ict, k12, p1_cloud, p2_cloud, k12, df, &
828  & spline_o%n1, spline_o%x1pkg(1,1), &
829  & spline_o%n2, spline_o%x2pkg(1,1), &
830  & spline_o%fspl(1,1,1), spline_o%n1, &
831  & iwarn, ifail)
832 
833  else if(spline_o%isHermite == 0) then
834 
835  call vecbicub(ict, k12, p1_cloud, p2_cloud, k12, df, &
836  & spline_o%n1,spline_o%x1pkg(1,1), &
837  & spline_o%n2,spline_o%x2pkg(1,1), &
838  & spline_o%fspl(1,1,1), spline_o%n1, &
839  & iwarn, ifail)
840 
841  else
842 
843  call vecherm2(ict, k12, p1_cloud, p2_cloud, k12, df, &
844  & spline_o%n1, spline_o%x1pkg(1,1), &
845  & spline_o%n2, spline_o%x2pkg(1,1), &
846  & spline_o%fspl(1,1,1), spline_o%n1, &
847  & iwarn, ifail)
848 
849  endif
850 
851  deallocate(p1_cloud, p2_cloud, stat=ifail)
852  if(ifail /= 0) then
853  ier = 31
854  return
855  endif
856 
857  if(ifail /= 0) ier = 95
858 
859 end subroutine ezspline_gradient2_array_r4
860 
861 
862 !!!
863 !!! 3-D
864 !!!
865 
866 subroutine ezspline_gradient3_r4(spline_o, &
867  p1, p2, p3, df, ier)
868  use ezspline_obj
869  implicit none
870  type(EZspline3_r4) spline_o
871  real(ezspline_r4), intent(in) :: p1, p2, p3
872  real(ezspline_r4), intent(out) :: df(3)
873  !
874  integer, intent(out) :: ier
875  integer ifail
876  integer, parameter :: ict(10) = (/0, 1, 1, 1, 0, 0, 0, 0, 0, 0/)
877 
878  ier = 0
879  ifail=0
880  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
881  ier = 94
882  return
883  endif
884 
885  if (spline_o%isHybrid == 1) then
886 
887  call evintrp3d(p1, p2, p3, &
888  & spline_o%x1(1), spline_o%n1, &
889  & spline_o%x2(1), spline_o%n2, &
890  & spline_o%x3(1), spline_o%n3, &
891  & spline_o%hspline, spline_o%fspl(1,1,1,1), &
892  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
893  & size(spline_o%fspl,3), size(spline_o%fspl,4), &
894  & ict, df, ifail)
895 
896  else if(spline_o%isLinear == 1) then
897 
898  call pc3ev(p1, p2, p3, &
899  & spline_o%x1(1), spline_o%n1, &
900  & spline_o%x2(1), spline_o%n2, &
901  & spline_o%x3(1), spline_o%n3, &
902  & spline_o%ilin1, spline_o%ilin2, spline_o%ilin3, &
903  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
904  & ict, df, ifail)
905 
906  else if(spline_o%isHermite == 0) then
907 
908  call evtricub(p1, p2, p3, &
909  & spline_o%x1(1), spline_o%n1, &
910  & spline_o%x2(1), spline_o%n2, &
911  & spline_o%x3(1), spline_o%n3, &
912  & spline_o%ilin1, spline_o%ilin2, spline_o%ilin3, &
913  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
914  & ict, df, ifail)
915 
916  else
917 
918  call herm3ev(p1, p2, p3, &
919  & spline_o%x1(1), spline_o%n1, &
920  & spline_o%x2(1), spline_o%n2, &
921  & spline_o%x3(1), spline_o%n3, &
922  & spline_o%ilin1, spline_o%ilin2, spline_o%ilin3, &
923  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
924  & ict, df, ifail)
925 
926  endif
927 
928  if(ifail/=0) ier = 95
929 
930 end subroutine ezspline_gradient3_r4
931 
932 subroutine ezspline_gradient3_cloud_r4(spline_o, k, &
933  p1, p2, p3, df, ier)
934  use ezspline_obj
935  implicit none
936  type(EZspline3_r4) spline_o
937  integer, intent(in) :: k
938  real(ezspline_r4), intent(in) :: p1(k), p2(k), p3(k)
939  real(ezspline_r4), intent(out) :: df(k, 3)
940  !
941  integer, intent(out) :: ier
942  integer ifail
943  integer, parameter :: ict(10) = (/0, 1, 1, 1, 0, 0, 0, 0, 0, 0/)
944  integer :: iwarn=0
945 
946  ier = 0
947  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
948  ier = 94
949  return
950  endif
951 
952  if (spline_o%isHybrid == 1) then
953 
954  call vecintrp3d(ict, k, p1, p2, p3, k, df, &
955  & spline_o%n1, spline_o%x1pkg(1,1), &
956  & spline_o%n2, spline_o%x2pkg(1,1), &
957  & spline_o%n3, spline_o%x3pkg(1,1), &
958  & spline_o%hspline, spline_o%fspl(1,1,1,1), &
959  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
960  & size(spline_o%fspl,3), size(spline_o%fspl,4), &
961  & iwarn, ifail)
962 
963  else if(spline_o%isLinear == 1) then
964 
965  call vecpc3(ict, k, p1, p2, p3, k, df, &
966  & spline_o%n1, spline_o%x1pkg(1,1), &
967  & spline_o%n2, spline_o%x2pkg(1,1), &
968  & spline_o%n3, spline_o%x3pkg(1,1), &
969  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
970  & iwarn,ifail)
971 
972  else if(spline_o%isHermite == 0) then
973 
974  call vectricub(ict, k, p1, p2, p3, k, df, &
975  & spline_o%n1,spline_o%x1pkg(1,1), &
976  & spline_o%n2,spline_o%x2pkg(1,1), &
977  & spline_o%n3,spline_o%x3pkg(1,1), &
978  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
979  & iwarn, ifail)
980 
981  else
982 
983  call vecherm3(ict, k, p1, p2, p3, k, df, &
984  & spline_o%n1, spline_o%x1pkg(1,1), &
985  & spline_o%n2, spline_o%x2pkg(1,1), &
986  & spline_o%n3, spline_o%x3pkg(1,1), &
987  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
988  & iwarn,ifail)
989 
990  endif
991 
992  if(ifail /= 0) ier = 95
993 
994 
995 end subroutine ezspline_gradient3_cloud_r4
996 
997 
998 subroutine ezspline_gradient3_array_r4(spline_o, k1, k2, k3, &
999  p1, p2, p3, df, ier)
1000  use ezspline_obj
1001  implicit none
1002  type(EZspline3_r4) spline_o
1003  integer, intent(in) :: k1, k2, k3
1004  real(ezspline_r4), intent(in) :: p1(k1), p2(k2), p3(k3)
1005  real(ezspline_r4), intent(out) :: df(k1, k2, k3, 3)
1006  !
1007  integer, intent(out) :: ier
1008  integer ifail
1009  integer, parameter :: ict(10) = (/0, 1, 1, 1, 0, 0, 0, 0, 0, 0/)
1010  integer:: iwarn=0
1011  real(ezspline_r4), dimension(:), allocatable :: p1_cloud, p2_cloud, p3_cloud
1012  integer k123
1013 
1014  ier = 0
1015  ifail=0
1016  if( .not.ezspline_allocated(spline_o) .or. spline_o%isReady /= 1) then
1017  ier = 94
1018  return
1019  endif
1020 
1021  k123 = k1*k2*k3
1022  allocate(p1_cloud(k123), p2_cloud(k123), p3_cloud(k123), stat=ifail)
1023  if(ifail /= 0) then
1024  ier = 30
1025  return
1026  endif
1027 
1028  p1_cloud = reshape(source=spread( &
1029  & source=spread(source=p1, dim=2, ncopies=k2), &
1030  & dim=3, ncopies=k3), shape=(/k123/))
1031  p2_cloud = reshape(source=spread( &
1032  & source=spread(source=p2, dim=1, ncopies=k1), &
1033  & dim=3, ncopies=k3), shape=(/k123/))
1034  p3_cloud = reshape(source=spread( &
1035  & source=spread(source=p3, dim=1, ncopies=k1), &
1036  & dim=2, ncopies=k2), shape=(/k123/))
1037 
1038  if (spline_o%isHybrid == 1) then
1039 
1040  call vecintrp3d(ict, k123, p1_cloud, p2_cloud, p3_cloud, k123, df, &
1041  & spline_o%n1, spline_o%x1pkg(1,1), &
1042  & spline_o%n2, spline_o%x2pkg(1,1), &
1043  & spline_o%n3, spline_o%x3pkg(1,1), &
1044  & spline_o%hspline, spline_o%fspl(1,1,1,1), &
1045  & size(spline_o%fspl,1), size(spline_o%fspl,2), &
1046  & size(spline_o%fspl,3), size(spline_o%fspl,4), &
1047  & iwarn, ifail)
1048 
1049  else if(spline_o%isLinear == 1) then
1050 
1051  call vecpc3(ict, k123, p1_cloud, p2_cloud, p3_cloud, k123, df, &
1052  & spline_o%n1, spline_o%x1pkg(1,1), &
1053  & spline_o%n2, spline_o%x2pkg(1,1), &
1054  & spline_o%n3, spline_o%x3pkg(1,1), &
1055  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
1056  & iwarn,ifail)
1057 
1058  else if(spline_o%isHermite == 0) then
1059 
1060  call vectricub(ict, k123, p1_cloud, p2_cloud, p3_cloud, k123, df, &
1061  & spline_o%n1,spline_o%x1pkg(1,1), &
1062  & spline_o%n2,spline_o%x2pkg(1,1), &
1063  & spline_o%n3,spline_o%x3pkg(1,1), &
1064  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
1065  & iwarn, ifail)
1066 
1067  else
1068 
1069  call vecherm3(ict, k123, p1_cloud, p2_cloud, p3_cloud, k123, df, &
1070  & spline_o%n1, spline_o%x1pkg(1,1), &
1071  & spline_o%n2, spline_o%x2pkg(1,1), &
1072  & spline_o%n3, spline_o%x3pkg(1,1), &
1073  & spline_o%fspl(1,1,1,1), spline_o%n1, spline_o%n2, &
1074  & iwarn,ifail)
1075 
1076  endif
1077 
1078  if(ifail /= 0) ier = 95
1079 
1080  deallocate(p1_cloud, p2_cloud, p3_cloud, stat=ifail)
1081  if(ifail /= 0) then
1082  ier = 31
1083  return
1084  endif
1085 
1086 end subroutine ezspline_gradient3_array_r4
ezspline_obj::ezspline_allocated
Definition: ezspline_obj.f90:19