V3FIT
ezhybrid_init.f90
1 !/////
2 ! R8 !
3 !/////
4 
5 subroutine ezhybrid_init2_r8(spline_o, n1, n2, hspline, ier, &
6  BCS1, BCS2)
7  use ezspline_obj
8  implicit none
9  type(EZspline2_r8) spline_o
10  integer, intent(in) :: n1, n2
11  integer, intent(in) :: hspline(2)
12  ! ier:
13  ! 0=ok
14  ! 1=allocation error
15  ! 99=something strange happened in EZhybrid_init
16  integer, intent(out) :: ier
17 
18  integer, intent(in), OPTIONAL :: BCS1(2), BCS2(2)
19 
20  integer i, iok, icoeff, idim1, idim2
21  logical :: spline_present, hermite_present
22 
23  ier = 0
24 
25  if(ezspline_allocated(spline_o)) then
26  ier = 100 ! allocated already
27  return
28  else
29  call ezspline_preinit(spline_o)
30  endif
31 
32  spline_o%n1 = n1
33  spline_o%n2 = n2
34 
35  spline_o%klookup1 = 3 ! default lookup algorithm selection
36  spline_o%klookup2 = 3 ! default lookup algorithm selection
37 
38  spline_o%isHermite = 0
39  spline_o%isLinear = 0
40  spline_o%isHybrid = 1
41  spline_o%hspline = hspline
42 
43  spline_present=.false.
44  hermite_present=.false.
45  icoeff=1
46  do i=1,2
47  if(hspline(i).lt.-1) then
48  ier=54
49  else if(hspline(i).gt.2) then
50  ier=54
51  else if(hspline(i).eq.1) then
52  hermite_present=.true.
53  icoeff=icoeff*2
54  else if(hspline(i).eq.2) then
55  spline_present=.true.
56  icoeff=icoeff*2
57  endif
58  enddo
59  if((ier.eq.0).and.spline_present.and.hermite_present) ier=55
60  if(ier.ne.0) return
61 
62  idim1 = n1
63  if(hspline(1).eq.-1) idim1 = n1 - 1
64  idim2 = n2
65  if(hspline(2).eq.-1) idim2 = n2 - 1
66 
67  iok = 0
68  allocate(spline_o%x1(n1), stat=iok)
69  if(iok /= 0) ier = 1
70  allocate(spline_o%x2(n2), stat=iok)
71  if(iok /= 0) ier = 1
72  allocate(spline_o%fspl(icoeff,idim1,idim2), stat=iok)
73  if(iok /= 0) ier = 1
74  allocate(spline_o%bcval1min(idim2), stat=iok)
75  if(iok /= 0) ier = 1
76  allocate(spline_o%bcval1max(idim2), stat=iok)
77  if(iok /= 0) ier = 1
78  allocate(spline_o%bcval2min(idim1), stat=iok)
79  if(iok /= 0) ier = 1
80  allocate(spline_o%bcval2max(idim1), stat=iok)
81  if(iok /= 0) ier = 1
82  allocate(spline_o%x1pkg(n1, 4), stat=iok)
83  if(iok /= 0) ier = 1
84  allocate(spline_o%x2pkg(n2, 4), stat=iok)
85  if(iok /= 0) ier = 1
86 
87  do i = 1, 2
88 
89  spline_o%bcval1min(1:idim2) = 0.0_ezspline_r8
90  spline_o%bcval1max(1:idim2) = 0.0_ezspline_r8
91  spline_o%ibctype1(i) = 0
92  if(present(bcs1)) then
93  select case(bcs1(i))
94  case (-1)
95  spline_o%ibctype1(i) = -1
96  case (0)
97  spline_o%ibctype1(i) = 0
98  case (1)
99  if(hspline(1).le.0) then
100  spline_o%ibctype1(i) = 0
101  ier = 56
102  else
103  spline_o%ibctype1(i) = 1
104  endif
105  case (2)
106  if(hspline(1).le.1) then
107  spline_o%ibctype1(i) = 0
108  ier = 56
109  else
110  spline_o%ibctype1(i) = 2
111  endif
112  case default
113  ier = 2
114  spline_o%ibctype1(i) = 0
115  end select
116  endif
117 
118  spline_o%bcval2min(1:idim1) = 0.0_ezspline_r8
119  spline_o%bcval2max(1:idim1) = 0.0_ezspline_r8
120  spline_o%ibctype2(i) = 0
121  if(present(bcs2)) then
122  select case(bcs2(i))
123  case (-1)
124  spline_o%ibctype2(i) = -1
125  case (0)
126  spline_o%ibctype2(i) = 0
127  case (1)
128  if(hspline(2).le.0) then
129  spline_o%ibctype2(i) = 0
130  ier = 56
131  else
132  spline_o%ibctype2(i) = 1
133  endif
134  case (2)
135  if(hspline(2).le.1) then
136  spline_o%ibctype2(i) = 0
137  ier = 56
138  else
139  spline_o%ibctype2(i) = 2
140  endif
141  case default
142  ier = 2
143  spline_o%ibctype2(i) = 0
144  end select
145  endif
146 
147  enddo
148 
149  if(spline_o%ibctype1(1)==-1 .OR. spline_o%ibctype1(2)==-1) then
150  spline_o%ibctype1(1)=-1
151  spline_o%ibctype1(2)=-1
152  endif
153  if(spline_o%ibctype2(1)==-1 .OR. spline_o%ibctype2(2)==-1) then
154  spline_o%ibctype2(1)=-1
155  spline_o%ibctype2(2)=-1
156  endif
157 
158  ! default grid
159  spline_o%x1min = 0.0_ezspline_r8
160  spline_o%x1max = 1.0_ezspline_r8
161  if(spline_o%ibctype1(1).eq.-1) spline_o%x1max = ezspline_twopi_r8
162 
163  spline_o%x1 = spline_o%x1min + (spline_o%x1max - spline_o%x1min)* &
164  & (/ (real(i-1,ezspline_r8)/real(spline_o%n1-1, ezspline_r8), i=1,spline_o%n1 ) /)
165  spline_o%x2min = 0.0_ezspline_r8
166  spline_o%x2max = 1.0_ezspline_r8
167  if(spline_o%ibctype2(1).eq.-1) spline_o%x2max = ezspline_twopi_r8
168  spline_o%x2 = spline_o%x2min + (spline_o%x2max - spline_o%x2min)* &
169  & (/ (real(i-1,ezspline_r8)/real(spline_o%n2-1, ezspline_r8), i=1,spline_o%n2 ) /)
170 
171  spline_o%isReady = 0
172 
173  return
174 end subroutine ezhybrid_init2_r8
175 
176 
177 
178 
179 subroutine ezhybrid_init3_r8(spline_o, n1, n2, n3, hspline, ier, &
180  BCS1, BCS2, BCS3)
181  use ezspline_obj
182  implicit none
183  type(EZspline3_r8) spline_o
184  integer, intent(in) :: n1, n2, n3
185  integer, intent(in) :: hspline(3)
186  ! ier:
187  ! 0=ok
188  ! 1=allocation error
189  ! 99=something strange happened in EZhybrid_init
190  integer, intent(out) :: ier
191 
192  integer, intent(in), OPTIONAL :: BCS1(2), BCS2(2), BCS3(2)
193 
194  integer i, iok, icoeff, idim1, idim2, idim3
195  logical :: spline_present, hermite_present
196 
197  ier = 0
198 
199  if(ezspline_allocated(spline_o)) then
200  ier = 100 ! allocated already
201  return
202  else
203  call ezspline_preinit(spline_o)
204  endif
205 
206  spline_o%n1 = n1
207  spline_o%n2 = n2
208  spline_o%n3 = n3
209 
210  spline_o%klookup1 = 3 ! default lookup algorithm selection
211  spline_o%klookup2 = 3 ! default lookup algorithm selection
212  spline_o%klookup3 = 3 ! default lookup algorithm selection
213 
214  spline_o%isHermite = 0
215  spline_o%isLinear = 0
216  spline_o%isHybrid = 1
217  spline_o%hspline = hspline
218 
219  spline_present=.false.
220  hermite_present=.false.
221  icoeff=1
222  do i=1,3
223  if(hspline(i).lt.-1) then
224  ier=54
225  else if(hspline(i).gt.2) then
226  ier=54
227  else if(hspline(i).eq.1) then
228  hermite_present=.true.
229  icoeff=icoeff*2
230  else if(hspline(i).eq.2) then
231  spline_present=.true.
232  icoeff=icoeff*2
233  endif
234  enddo
235  if((ier.eq.0).and.spline_present.and.hermite_present) ier=55
236  if(ier.ne.0) return
237 
238  idim1 = n1
239  if(hspline(1).eq.-1) idim1 = n1 - 1
240  idim2 = n2
241  if(hspline(2).eq.-1) idim2 = n2 - 1
242  idim3 = n3
243  if(hspline(3).eq.-1) idim3 = n3 - 1
244 
245  iok = 0
246  allocate(spline_o%x1(n1), stat=iok)
247  if(iok /= 0) ier = 1
248  allocate(spline_o%x2(n2), stat=iok)
249  if(iok /= 0) ier = 1
250  allocate(spline_o%x3(n3), stat=iok)
251  if(iok /= 0) ier = 1
252  allocate(spline_o%fspl(icoeff,idim1,idim2,idim3), stat=iok)
253  if(iok /= 0) ier = 1
254  allocate(spline_o%bcval1min(idim2, idim3), stat=iok)
255  if(iok /= 0) ier = 1
256  allocate(spline_o%bcval1max(idim2, idim3), stat=iok)
257  if(iok /= 0) ier = 1
258  allocate(spline_o%bcval2min(idim1, idim3), stat=iok)
259  if(iok /= 0) ier = 1
260  allocate(spline_o%bcval2max(idim1, idim3), stat=iok)
261  if(iok /= 0) ier = 1
262  allocate(spline_o%bcval3min(idim1, idim2), stat=iok)
263  if(iok /= 0) ier = 1
264  allocate(spline_o%bcval3max(idim1, idim2), stat=iok)
265  if(iok /= 0) ier = 1
266  allocate(spline_o%x1pkg(n1, 4), stat=iok)
267  if(iok /= 0) ier = 1
268  allocate(spline_o%x2pkg(n2, 4), stat=iok)
269  if(iok /= 0) ier = 1
270  allocate(spline_o%x3pkg(n3, 4), stat=iok)
271  if(iok /= 0) ier = 1
272 
273  do i = 1, 2
274 
275  spline_o%bcval1min(1:idim2, 1:idim3) = 0.0_ezspline_r8
276  spline_o%bcval1max(1:idim2, 1:idim3) = 0.0_ezspline_r8
277  spline_o%ibctype1(i) = 0
278  if(present(bcs1)) then
279  select case(bcs1(i))
280  case (-1)
281  spline_o%ibctype1(i) = -1
282  case (0)
283  spline_o%ibctype1(i) = 0
284  case (1)
285  if(hspline(1).le.0) then
286  spline_o%ibctype1(i) = 0
287  ier = 56
288  else
289  spline_o%ibctype1(i) = 1
290  endif
291  case (2)
292  if(hspline(1).le.1) then
293  spline_o%ibctype1(i) = 0
294  ier = 56
295  else
296  spline_o%ibctype1(i) = 2
297  endif
298  case default
299  ier = 2
300  spline_o%ibctype1(i) = 0
301  end select
302  endif
303 
304  spline_o%bcval2min(1:idim1, 1:idim3) = 0.0_ezspline_r8
305  spline_o%bcval2max(1:idim1, 1:idim3) = 0.0_ezspline_r8
306  spline_o%ibctype2(i) = 0
307  if(present(bcs2)) then
308  select case(bcs2(i))
309  case (-1)
310  spline_o%ibctype2(i) = -1
311  case (0)
312  spline_o%ibctype2(i) = 0
313  case (1)
314  if(hspline(2).le.0) then
315  spline_o%ibctype2(i) = 0
316  ier = 56
317  else
318  spline_o%ibctype2(i) = 1
319  endif
320  case (2)
321  if(hspline(2).le.1) then
322  spline_o%ibctype2(i) = 0
323  ier = 56
324  else
325  spline_o%ibctype2(i) = 2
326  endif
327  case default
328  ier = 2
329  spline_o%ibctype2(i) = 0
330  end select
331  endif
332 
333  spline_o%bcval3min(1:idim1, 1:idim2) = 0.0_ezspline_r8
334  spline_o%bcval3max(1:idim1, 1:idim2) = 0.0_ezspline_r8
335  spline_o%ibctype3(i) = 0
336  if(present(bcs3)) then
337  select case(bcs3(i))
338  case (-1)
339  spline_o%ibctype3(i) = -1
340  case (0)
341  spline_o%ibctype3(i) = 0
342  case (1)
343  if(hspline(3).le.0) then
344  spline_o%ibctype3(i) = 0
345  ier = 2
346  else
347  spline_o%ibctype3(i) = 1
348  endif
349  case (2)
350  if(hspline(3).le.1) then
351  spline_o%ibctype3(i) = 0
352  ier = 2
353  else
354  spline_o%ibctype3(i) = 2
355  endif
356  case default
357  ier = 2
358  spline_o%ibctype3(i) = 0
359  end select
360  endif
361 
362  enddo
363 
364  if(spline_o%ibctype1(1)==-1 .OR. spline_o%ibctype1(2)==-1) then
365  spline_o%ibctype1(1)=-1
366  spline_o%ibctype1(2)=-1
367  endif
368  if(spline_o%ibctype2(1)==-1 .OR. spline_o%ibctype2(2)==-1) then
369  spline_o%ibctype2(1)=-1
370  spline_o%ibctype2(2)=-1
371  endif
372  if(spline_o%ibctype3(1)==-1 .OR. spline_o%ibctype3(2)==-1) then
373  spline_o%ibctype3(1)=-1
374  spline_o%ibctype3(2)=-1
375  endif
376 
377  ! default grid
378  spline_o%x1min = 0.0_ezspline_r8
379  spline_o%x1max = 1.0_ezspline_r8
380  if(spline_o%ibctype1(1).eq.-1) spline_o%x1max = ezspline_twopi_r8
381  spline_o%x1 = spline_o%x1min + (spline_o%x1max - spline_o%x1min)* &
382  & (/ (real(i-1,ezspline_r8)/real(spline_o%n1-1, ezspline_r8), i=1,spline_o%n1 ) /)
383 
384  spline_o%x2min = 0.0_ezspline_r8
385  spline_o%x2max = 1.0_ezspline_r8
386  if(spline_o%ibctype2(1).eq.-1) spline_o%x2max = ezspline_twopi_r8
387  spline_o%x2 = spline_o%x2min + (spline_o%x2max - spline_o%x2min)* &
388  & (/ (real(i-1,ezspline_r8)/real(spline_o%n2-1, ezspline_r8), i=1,spline_o%n2 ) /)
389 
390  spline_o%x3min = 0.0_ezspline_r8
391  spline_o%x3max = 1.0_ezspline_r8
392  if(spline_o%ibctype3(1).eq.-1) spline_o%x3max = ezspline_twopi_r8
393  spline_o%x3 = spline_o%x3min + (spline_o%x3max - spline_o%x3min)* &
394  & (/ (real(i-1,ezspline_r8)/real(spline_o%n3-1, ezspline_r8), i=1,spline_o%n3 ) /)
395 
396  spline_o%isReady = 0
397 
398  return
399 end subroutine ezhybrid_init3_r8
400 
401 !/////
402 ! R4 !
403 !/////
404 
405 subroutine ezhybrid_init2_r4(spline_o, n1, n2, hspline, ier, &
406  BCS1, BCS2)
407  use ezspline_obj
408  implicit none
409  type(EZspline2_r4) spline_o
410  integer, intent(in) :: n1, n2
411  integer, intent(in) :: hspline(2)
412  ! ier:
413  ! 0=ok
414  ! 1=allocation error
415  ! 99=something strange happened in EZhybrid_init
416  integer, intent(out) :: ier
417 
418  integer, intent(in), OPTIONAL :: BCS1(2), BCS2(2)
419 
420  integer i, iok, icoeff, idim1, idim2
421  logical :: spline_present, hermite_present
422 
423  ier = 0
424 
425  if(ezspline_allocated(spline_o)) then
426  ier = 100 ! allocated already
427  return
428  else
429  call ezspline_preinit(spline_o)
430  endif
431 
432  spline_o%n1 = n1
433  spline_o%n2 = n2
434 
435  spline_o%klookup1 = 3 ! default lookup algorithm selection
436  spline_o%klookup2 = 3 ! default lookup algorithm selection
437 
438  spline_o%isHermite = 0
439  spline_o%isLinear = 0
440  spline_o%isHybrid = 1
441  spline_o%hspline = hspline
442 
443  spline_present=.false.
444  hermite_present=.false.
445  icoeff=1
446  do i=1,2
447  if(hspline(i).lt.-1) then
448  ier=54
449  else if(hspline(i).gt.2) then
450  ier=54
451  else if(hspline(i).eq.1) then
452  hermite_present=.true.
453  icoeff=icoeff*2
454  else if(hspline(i).eq.2) then
455  spline_present=.true.
456  icoeff=icoeff*2
457  endif
458  enddo
459  if((ier.eq.0).and.spline_present.and.hermite_present) ier=55
460  if(ier.ne.0) return
461 
462  idim1 = n1
463  if(hspline(1).eq.-1) idim1 = n1 - 1
464  idim2 = n2
465  if(hspline(2).eq.-1) idim2 = n2 - 1
466 
467  iok = 0
468  allocate(spline_o%x1(n1), stat=iok)
469  if(iok /= 0) ier = 1
470  allocate(spline_o%x2(n2), stat=iok)
471  if(iok /= 0) ier = 1
472  allocate(spline_o%fspl(icoeff,idim1,idim2), stat=iok)
473  if(iok /= 0) ier = 1
474  allocate(spline_o%bcval1min(idim2), stat=iok)
475  if(iok /= 0) ier = 1
476  allocate(spline_o%bcval1max(idim2), stat=iok)
477  if(iok /= 0) ier = 1
478  allocate(spline_o%bcval2min(idim1), stat=iok)
479  if(iok /= 0) ier = 1
480  allocate(spline_o%bcval2max(idim1), stat=iok)
481  if(iok /= 0) ier = 1
482  allocate(spline_o%x1pkg(n1, 4), stat=iok)
483  if(iok /= 0) ier = 1
484  allocate(spline_o%x2pkg(n2, 4), stat=iok)
485  if(iok /= 0) ier = 1
486 
487  do i = 1, 2
488 
489  spline_o%bcval1min(1:idim2) = 0.0
490  spline_o%bcval1max(1:idim2) = 0.0
491  spline_o%ibctype1(i) = 0
492  if(present(bcs1)) then
493  select case(bcs1(i))
494  case (-1)
495  spline_o%ibctype1(i) = -1
496  case (0)
497  spline_o%ibctype1(i) = 0
498  case (1)
499  if(hspline(1).le.0) then
500  spline_o%ibctype1(i) = 0
501  ier = 56
502  else
503  spline_o%ibctype1(i) = 1
504  endif
505  case (2)
506  if(hspline(1).le.1) then
507  spline_o%ibctype1(i) = 0
508  ier = 56
509  else
510  spline_o%ibctype1(i) = 2
511  endif
512  case default
513  ier = 2
514  spline_o%ibctype1(i) = 0
515  end select
516  endif
517 
518  spline_o%bcval2min(1:idim1) = 0.0
519  spline_o%bcval2max(1:idim1) = 0.0
520  spline_o%ibctype2(i) = 0
521  if(present(bcs2)) then
522  select case(bcs2(i))
523  case (-1)
524  spline_o%ibctype2(i) = -1
525  case (0)
526  spline_o%ibctype2(i) = 0
527  case (1)
528  if(hspline(2).le.0) then
529  spline_o%ibctype2(i) = 0
530  ier = 56
531  else
532  spline_o%ibctype2(i) = 1
533  endif
534  case (2)
535  if(hspline(2).le.1) then
536  spline_o%ibctype2(i) = 0
537  ier = 56
538  else
539  spline_o%ibctype2(i) = 2
540  endif
541  case default
542  ier = 2
543  spline_o%ibctype2(i) = 0
544  end select
545  endif
546 
547  enddo
548 
549  if(spline_o%ibctype1(1)==-1 .OR. spline_o%ibctype1(2)==-1) then
550  spline_o%ibctype1(1)=-1
551  spline_o%ibctype1(2)=-1
552  endif
553  if(spline_o%ibctype2(1)==-1 .OR. spline_o%ibctype2(2)==-1) then
554  spline_o%ibctype2(1)=-1
555  spline_o%ibctype2(2)=-1
556  endif
557 
558  ! default grid
559  spline_o%x1min = 0.0_ezspline_r4
560  spline_o%x1max = 1.0_ezspline_r4
561  if(spline_o%ibctype1(1).eq.-1) spline_o%x1max = ezspline_twopi_r4
562  spline_o%x1 = spline_o%x1min + (spline_o%x1max - spline_o%x1min)* &
563  & (/ (real(i-1,ezspline_r4)/real(spline_o%n1-1, ezspline_r4), i=1,spline_o%n1 ) /)
564 
565  spline_o%x2min = 0.0_ezspline_r4
566  spline_o%x2max = 1.0_ezspline_r4
567  if(spline_o%ibctype2(1).eq.-1) spline_o%x2max = ezspline_twopi_r4
568  spline_o%x2 = spline_o%x2min + (spline_o%x2max - spline_o%x2min)* &
569  & (/ (real(i-1,ezspline_r4)/real(spline_o%n2-1, ezspline_r4), i=1,spline_o%n2 ) /)
570 
571 
572  spline_o%isReady = 0
573 
574  return
575 end subroutine ezhybrid_init2_r4
576 
577 
578 
579 
580 subroutine ezhybrid_init3_r4(spline_o, n1, n2, n3, hspline, ier, &
581  BCS1, BCS2, BCS3)
582  use ezspline_obj
583  implicit none
584  type(EZspline3_r4) spline_o
585  integer, intent(in) :: n1, n2, n3
586  integer, intent(in) :: hspline(3)
587  ! ier:
588  ! 0=ok
589  ! 1=allocation error
590  ! 99=something strange happened in EZhybrid_init
591  integer, intent(out) :: ier
592 
593  integer, intent(in), OPTIONAL :: BCS1(2), BCS2(2), BCS3(2)
594 
595  integer i, iok, icoeff, idim1, idim2, idim3
596  logical :: spline_present, hermite_present
597 
598  ier = 0
599 
600  if(ezspline_allocated(spline_o)) then
601  ier = 100 ! allocated already
602  return
603  else
604  call ezspline_preinit(spline_o)
605  endif
606 
607  spline_o%n1 = n1
608  spline_o%n2 = n2
609  spline_o%n3 = n3
610 
611  spline_o%klookup1 = 3 ! default lookup algorithm selection
612  spline_o%klookup2 = 3 ! default lookup algorithm selection
613  spline_o%klookup3 = 3 ! default lookup algorithm selection
614 
615  spline_o%isHermite = 0
616  spline_o%isLinear = 0
617  spline_o%isHybrid = 1
618  spline_o%hspline = hspline
619 
620  spline_present=.false.
621  hermite_present=.false.
622  icoeff=1
623  do i=1,3
624  if(hspline(i).lt.-1) then
625  ier=54
626  else if(hspline(i).gt.2) then
627  ier=54
628  else if(hspline(i).eq.1) then
629  hermite_present=.true.
630  icoeff=icoeff*2
631  else if(hspline(i).eq.2) then
632  spline_present=.true.
633  icoeff=icoeff*2
634  endif
635  enddo
636  if((ier.eq.0).and.spline_present.and.hermite_present) ier=55
637  if(ier.ne.0) return
638 
639  idim1 = n1
640  if(hspline(1).eq.-1) idim1 = n1 - 1
641  idim2 = n2
642  if(hspline(2).eq.-1) idim2 = n2 - 1
643  idim3 = n3
644  if(hspline(3).eq.-1) idim3 = n3 - 1
645 
646  iok = 0
647  allocate(spline_o%x1(n1), stat=iok)
648  if(iok /= 0) ier = 1
649  allocate(spline_o%x2(n2), stat=iok)
650  if(iok /= 0) ier = 1
651  allocate(spline_o%x3(n3), stat=iok)
652  if(iok /= 0) ier = 1
653  allocate(spline_o%fspl(icoeff,idim1,idim2,idim3), stat=iok)
654  if(iok /= 0) ier = 1
655  allocate(spline_o%bcval1min(idim2, idim3), stat=iok)
656  if(iok /= 0) ier = 1
657  allocate(spline_o%bcval1max(idim2, idim3), stat=iok)
658  if(iok /= 0) ier = 1
659  allocate(spline_o%bcval2min(idim1, idim3), stat=iok)
660  if(iok /= 0) ier = 1
661  allocate(spline_o%bcval2max(idim1, idim3), stat=iok)
662  if(iok /= 0) ier = 1
663  allocate(spline_o%bcval3min(idim1, idim2), stat=iok)
664  if(iok /= 0) ier = 1
665  allocate(spline_o%bcval3max(idim1, idim2), stat=iok)
666  if(iok /= 0) ier = 1
667  allocate(spline_o%x1pkg(n1, 4), stat=iok)
668  if(iok /= 0) ier = 1
669  allocate(spline_o%x2pkg(n2, 4), stat=iok)
670  if(iok /= 0) ier = 1
671  allocate(spline_o%x3pkg(n3, 4), stat=iok)
672  if(iok /= 0) ier = 1
673 
674  do i = 1, 2
675 
676  spline_o%bcval1min(1:idim2, 1:idim3) = 0.0
677  spline_o%bcval1max(1:idim2, 1:idim3) = 0.0
678  spline_o%ibctype1(i) = 0
679  if(present(bcs1)) then
680  select case(bcs1(i))
681  case (-1)
682  spline_o%ibctype1(i) = -1
683  case (0)
684  spline_o%ibctype1(i) = 0
685  case (1)
686  if(hspline(1).le.0) then
687  spline_o%ibctype1(i) = 0
688  ier = 56
689  else
690  spline_o%ibctype1(i) = 1
691  endif
692  case (2)
693  if(hspline(1).le.1) then
694  spline_o%ibctype1(i) = 0
695  ier = 56
696  else
697  spline_o%ibctype1(i) = 2
698  endif
699  case default
700  ier = 2
701  spline_o%ibctype1(i) = 0
702  end select
703  endif
704 
705  spline_o%bcval2min(1:idim1, 1:idim3) = 0.0
706  spline_o%bcval2max(1:idim1, 1:idim3) = 0.0
707  spline_o%ibctype2(i) = 0
708  if(present(bcs2)) then
709  select case(bcs2(i))
710  case (-1)
711  spline_o%ibctype2(i) = -1
712  case (0)
713  spline_o%ibctype2(i) = 0
714  case (1)
715  if(hspline(2).le.0) then
716  spline_o%ibctype2(i) = 0
717  ier = 56
718  else
719  spline_o%ibctype2(i) = 1
720  endif
721  case (2)
722  if(hspline(2).le.1) then
723  spline_o%ibctype2(i) = 0
724  ier = 56
725  else
726  spline_o%ibctype2(i) = 2
727  endif
728  case default
729  ier = 2
730  spline_o%ibctype2(i) = 0
731  end select
732  endif
733 
734  spline_o%bcval3min(1:idim1, 1:idim2) = 0.0
735  spline_o%bcval3max(1:idim1, 1:idim2) = 0.0
736  spline_o%ibctype3(i) = 0
737  if(present(bcs3)) then
738  select case(bcs3(i))
739  case (-1)
740  spline_o%ibctype3(i) = -1
741  case (0)
742  spline_o%ibctype3(i) = 0
743  case (1)
744  if(hspline(3).le.0) then
745  spline_o%ibctype3(i) = 0
746  ier = 2
747  else
748  spline_o%ibctype3(i) = 1
749  endif
750  case (2)
751  if(hspline(3).le.1) then
752  spline_o%ibctype3(i) = 0
753  ier = 2
754  else
755  spline_o%ibctype3(i) = 2
756  endif
757  case default
758  ier = 2
759  spline_o%ibctype3(i) = 0
760  end select
761  endif
762 
763  enddo
764 
765  if(spline_o%ibctype1(1)==-1 .OR. spline_o%ibctype1(2)==-1) then
766  spline_o%ibctype1(1)=-1
767  spline_o%ibctype1(2)=-1
768  endif
769  if(spline_o%ibctype2(1)==-1 .OR. spline_o%ibctype2(2)==-1) then
770  spline_o%ibctype2(1)=-1
771  spline_o%ibctype2(2)=-1
772  endif
773  if(spline_o%ibctype3(1)==-1 .OR. spline_o%ibctype3(2)==-1) then
774  spline_o%ibctype3(1)=-1
775  spline_o%ibctype3(2)=-1
776  endif
777 
778  ! default grid
779  spline_o%x1min = 0.0_ezspline_r4
780  spline_o%x1max = 1.0_ezspline_r4
781  if(spline_o%ibctype1(1).eq.-1) spline_o%x1max = ezspline_twopi_r4
782  spline_o%x1 = spline_o%x1min + (spline_o%x1max - spline_o%x1min)* &
783  & (/ (real(i-1,ezspline_r4)/real(spline_o%n1-1, ezspline_r4), i=1,spline_o%n1 ) /)
784 
785  spline_o%x2min = 0.0_ezspline_r4
786  spline_o%x2max = 1.0_ezspline_r4
787  if(spline_o%ibctype2(1).eq.-1) spline_o%x2max = ezspline_twopi_r4
788  spline_o%x2 = spline_o%x2min + (spline_o%x2max - spline_o%x2min)* &
789  & (/ (real(i-1,ezspline_r4)/real(spline_o%n2-1, ezspline_r4), i=1,spline_o%n2 ) /)
790 
791  spline_o%x3min = 0.0_ezspline_r4
792  spline_o%x3max = 1.0_ezspline_r4
793  if(spline_o%ibctype3(1).eq.-1) spline_o%x3max = ezspline_twopi_r4
794  spline_o%x3 = spline_o%x3min + (spline_o%x3max - spline_o%x3min)* &
795  & (/ (real(i-1,ezspline_r4)/real(spline_o%n3-1, ezspline_r4), i=1,spline_o%n3 ) /)
796 
797  spline_o%isReady = 0
798 
799  return
800 end subroutine ezhybrid_init3_r4
801 
ezspline_obj::ezspline_allocated
Definition: ezspline_obj.f90:19
ezspline_obj::ezspline_preinit
Definition: ezspline_obj.f90:3