ETS  \$Id: Doxyfile 2162 2020-02-26 14:16:09Z g2dpc $
 All Classes Files Functions Variables Pages
eq_test.F90
Go to the documentation of this file.
1 ! + + + + + + + + + + + + + + + + + + + + + + + + + + + +
9 ! + + + + + + + + + + + + + + + + + + + + + + + + + + + +
10 program eq_test
11 
12  use itm_constants
13  use euitm_routines
14  USE euitm_schemas
15  use euitm_routines
16 #ifdef GOT_E3EMEQ
18 #endif
19 #ifdef GOT_HELENA
21 #endif
22 #ifdef GOT_CHEASE
24 #endif
25 #ifdef GOT_SPIDER
27 #endif
29  use ets_version
30  use xml_file_reader
31  use write_structures
32 
33  implicit none
34 
35  integer, parameter :: nrho=100
36  integer :: nrho1
37 
38  TYPE (type_equilibrium), POINTER :: eq_old(:), eq_new(:), eq_helena(:), eq_e3astr(:), eq_chease(:), eq_spider(:)
39  TYPE (type_coreprof), POINTER :: coreprof(:)
40  TYPE (type_param) :: code_parameters
41  integer :: irho
42  real (R8) :: a0, r0, z0, b0, rho
43  integer eq_code
44  integer :: shot_no
45  integer idx
46  integer iargc
47  character*256 arg
48  integer case
49  integer, parameter :: no_of_cases=4
50  real (R8) :: case_p0(no_of_cases)
51  integer case_shot(no_of_cases)
52 
53  data case_p0 /5000.0_r8, 500.0_r8, 50.0_r8, 5.0_r8/
54  data case_shot /7, 8, 9, 10/
55 
56  interface
57  SUBROUTINE helena21itm (euitm_eq_in, euitm_eq_out)
58  USE euitm_schemas
59  IMPLICIT NONE
60  type (type_equilibrium), pointer :: euitm_eq_in(:)
61  type (type_equilibrium), pointer :: euitm_eq_out(:)
62  end SUBROUTINE helena21itm
63  end interface
64 
65  if(iargc().eq.0) then
66  case = 1
67  else if(iargc().eq.1) then
68  call getarg(1,arg)
69  read(arg,*) case
70  if(case.lt.1.or.case.gt.no_of_cases) then
71  write(*,*) 'Only cases 1 to ',no_of_cases,' have been coded'
72  stop 'Case index out of range'
73  endif
74  else
75  write(*,*) 'Too many arguments'
76  stop 'Too many arguments'
77  endif
78 
79  write(*,*) 'Doing case ', case
80 
81 ! call perfinit
82 ! call perfon ('whole')
83  allocate(eq_old(1))
84  allocate(eq_old(1)%profiles_1d%rho_tor(nrho))
85  allocate(eq_old(1)%profiles_1d%pressure(nrho))
86  allocate(eq_old(1)%profiles_1d%jparallel(nrho))
87  allocate(eq_old(1)%profiles_1d%jphi(nrho))
88 
89  allocate(coreprof(1))
90  allocate(coreprof(1)%rho_tor(nrho))
91  allocate(coreprof(1)%te%value(nrho))
92 
93 ! initialize plasma
94 
95  a0 = 0.1_r8
96  r0 = 10.0_r8
97  z0 = 0.0_r8
98  b0 = 2.0_r8
99 
100  eq_old(1)%eqgeometry%geom_axis%r = r0
101  eq_old(1)%eqgeometry%geom_axis%z = z0
102  eq_old(1)%eqgeometry%a_minor = a0
103  eq_old(1)%global_param%mag_axis%bphi = b0
104  eq_old(1)%global_param%mag_axis%position%r = r0
105  eq_old(1)%global_param%mag_axis%position%Z = z0
106  eq_old(1)%global_param%toroid_field%b0 = b0
107  eq_old(1)%global_param%toroid_field%r0 = r0
108 
109  do irho=1, nrho
110  rho = real((irho-1), kind=r8)/real((nrho-1), kind=r8)
111  eq_old(1)%profiles_1d%rho_tor(irho) = rho*a0
112  eq_old(1)%profiles_1d%pressure(irho) = (1.0_r8 - rho**2.0_r8) * case_p0(case)
113  eq_old(1)%profiles_1d%jparallel(irho) = (1.0_r8 - rho**2.0_r8) * 1e4_r8 / (a0**2 * itm_pi)
114  coreprof(1)%rho_tor(irho) = rho*a0
115  coreprof(1)%te%value(irho) = sqrt(1.0_r8 - rho**2.0_r8) * case_p0(case)/1e19_r8/itm_ev
116  enddo
117  eq_old(1)%profiles_1d%jphi=eq_old(1)%profiles_1d%jparallel
118 
119  shot_no = case_shot(case)
120 
121  call open_write_file(10, 'eq_test.in.cpo')
122  call write_cpo(eq_old(1), 'equilibrium')
123  call close_write_file
124  open(10,file='eq_test.in')
125  write(10,'(''# rho pressure jparallel jphi'')')
126  do irho=1,nrho
127  write(10,'(30f20.6)') &
128  eq_old(1)%profiles_1d%rho_tor(irho), &
129  eq_old(1)%profiles_1d%pressure(irho), &
130  eq_old(1)%profiles_1d%jparallel(irho), &
131  eq_old(1)%profiles_1d%jphi(irho)
132  enddo
133  close(10)
134 
135 ! call bdseq (twice for testing) and then write results to 'eq_test.out'
136 
137 ! call perfon ('bdseq1')
138  call bdseq_wrapper(eq_old, eq_new)
139 ! call perfoff
140  call euitm_deallocate(eq_new)
141 ! call perfon ('bdseq2')
142 ! call perfoff
143  call bdseq_wrapper(eq_old, eq_new)
144  write(*,*) 'CODE: BDSEQ'
145  call output_eq_info(eq_new(1))
146  write(*,*)
147 
148  call open_write_file(10, 'eq_test.out.cpo')
149  call write_cpo(eq_new(1), 'equilibrium')
150  call close_write_file
151  open(10,file='eq_test.out')
152  write(10,'(''# rho pressure jparallel jphi volume q rhovol psi phi vprime gm1 gm2 gm3 gm4 gm5 gm6 gm7 pprime ffprime F_dia'')')
153  do irho=1,size(eq_new(1)%profiles_1d%rho_vol)
154  write(10,'(30f20.6)') &
155  eq_new(1)%profiles_1d%rho_tor(irho), &
156  eq_new(1)%profiles_1d%pressure(irho), &
157  eq_new(1)%profiles_1d%jparallel(irho), &
158  eq_new(1)%profiles_1d%jphi(irho), &
159  eq_new(1)%profiles_1d%volume(irho), &
160  eq_new(1)%profiles_1d%q(irho), &
161  eq_new(1)%profiles_1d%rho_vol(irho), &
162  eq_new(1)%profiles_1d%psi(irho), &
163  eq_new(1)%profiles_1d%phi(irho), &
164  eq_new(1)%profiles_1d%vprime(irho), &
165  eq_new(1)%profiles_1d%gm1(irho), &
166  eq_new(1)%profiles_1d%gm2(irho), &
167  eq_new(1)%profiles_1d%gm3(irho), &
168  eq_new(1)%profiles_1d%gm4(irho), &
169  eq_new(1)%profiles_1d%gm5(irho), &
170  eq_new(1)%profiles_1d%gm6(irho), &
171  eq_new(1)%profiles_1d%gm7(irho), &
172  eq_new(1)%profiles_1d%pprime(irho), &
173  eq_new(1)%profiles_1d%ffprime(irho), &
174  eq_new(1)%profiles_1d%F_dia(irho)
175  enddo
176  close(10)
177 
178  eq_new(1)%eqconstraint%bvac_r%measured = &
179  eq_old(1)%global_param%mag_axis%bphi * &
180  eq_old(1)%global_param%mag_axis%position%r
181  call euitm_deallocate(eq_old)
182 
183 #ifdef UAL
184 ! for BDSEQ, run_no = 1
185 ! call perfon ('ualw1')
186  write(*,*) 'Calling euitm_create for BDSEQ'
187  call euitm_create('euitm',shot_no,1,0,0,idx)
188  write(*,*) eq_new(1)%time
189  eq_new(1)%time=0.0_r8
190  if(.not.associated(eq_new(1)%codeparam%codename)) then
191  allocate(eq_new(1)%codeparam%codename(1))
192  eq_new(1)%codeparam%codename(1)='BDSEQ'
193  endif
194  if(.not.associated(eq_new(1)%codeparam%codeversion)) then
195  allocate(eq_new(1)%codeparam%codeversion(1))
196  eq_new(1)%codeparam%codeversion(1)=version
197  endif
198  if(.not.associated(eq_new(1)%codeparam%parameters)) then
199  allocate(eq_new(1)%codeparam%parameters(1))
200  eq_new(1)%codeparam%parameters(1)='NONE'
201  endif
202  if(.not.associated(eq_new(1)%codeparam%output_diag)) then
203  allocate(eq_new(1)%codeparam%output_diag(1))
204  eq_new(1)%codeparam%output_diag(1)='NONE'
205  endif
206  eq_new(1)%codeparam%output_flag=0
207 ! write(*,*) 'Calling euitm_put_non_timed for BDSEQ'
208 ! call euitm_put_non_timed(idx,"equilibrium",EQ_NEW(1))
209 ! write(*,*) 'Calling euitm_put_slice for BDSEQ'
210 ! call euitm_put_slice(idx,"equilibrium",EQ_NEW(1))
211  call euitm_put(idx,"equilibrium",eq_new)
212  write(*,*) 'Calling euitm_close for BDSEQ'
213  call euitm_close(idx)
214  write(*,*) 'Finished with euitm UAL for BDSEQ'
215 ! call perfoff
216 #endif
217 
218  do eq_code = 1,5
219 
220  if(eq_code.eq.1) then
221 
222 #ifdef GOT_HELENA
223 
224 ! call helena and write results to 'eq_helena.out'
225 
226 !!! now in wrapper EQ_NEW(1)%profiles_1d%jphi=-EQ_NEW(1)%profiles_1d%jphi
227 ! call perfon ('helena')
228  call helena_wrapper(eq_new, eq_helena)
229 ! call perfoff
230 !!! now in wrapper EQ_NEW(1)%profiles_1d%jphi=-EQ_NEW(1)%profiles_1d%jphi
231  write(*,*) 'CODE: HELENA'
232  call output_eq_info(eq_helena(1))
233  write(*,*)
234 
235  call open_write_file(10, 'eq_helena.out.cpo')
236  call write_cpo(eq_helena(1), 'equilibrium')
237  call close_write_file
238  open(10,file='eq_helena.out')
239 ! write(10,'(''# rho pressure jparallel jphi volume q rhovol psi phi vprime gm1 gm2 gm3 gm4 gm5 gm6 gm7'')')
240  write(10,'(''# pressure jphi q psi volume rho_vol F_dia pprime ffprime area rho_tor phi vprime gm1 gm2 gm3 gm4 gm5 gm6 gm7'')')
241  do irho=1,size(eq_helena(1)%profiles_1d%psi)
242  write(10,'(30f20.6)') &
243  eq_helena(1)%profiles_1d%pressure(irho), &
244  eq_helena(1)%profiles_1d%jphi(irho), &
245  eq_helena(1)%profiles_1d%q(irho), &
246  eq_helena(1)%profiles_1d%psi(irho), &
247  eq_helena(1)%profiles_1d%volume(irho), &
248  eq_helena(1)%profiles_1d%rho_vol(irho), &
249  eq_helena(1)%profiles_1d%F_dia(irho), &
250  eq_helena(1)%profiles_1d%pprime(irho), &
251  eq_helena(1)%profiles_1d%ffprime(irho), &
252  eq_helena(1)%profiles_1d%area(irho), &
253  eq_helena(1)%profiles_1d%rho_tor(irho), &
254  eq_helena(1)%profiles_1d%phi(irho), &
255  eq_helena(1)%profiles_1d%vprime(irho), &
256  eq_helena(1)%profiles_1d%gm1(irho), &
257  eq_helena(1)%profiles_1d%gm2(irho), &
258  eq_helena(1)%profiles_1d%gm3(irho), &
259  eq_helena(1)%profiles_1d%gm4(irho), &
260  eq_helena(1)%profiles_1d%gm5(irho), &
261  eq_helena(1)%profiles_1d%gm6(irho), &
262  eq_helena(1)%profiles_1d%gm7(irho)
263 ! EQ_HELENA(1)%profiles_1d%jparallel(irho), &
264  enddo
265  close(10)
266 
267 #ifdef UAL
268 ! for HELENA, run_no = 3
269 ! call perfon ('ualw2')
270  write(*,*) 'Calling euitm_create for HELENA'
271  call euitm_create('euitm',shot_no,3,0,0,idx)
272  write(*,*) eq_helena(1)%time
273  eq_helena(1)%time=0.0_r8
274  if(.not.associated(eq_helena(1)%codeparam%codename)) then
275  allocate(eq_helena(1)%codeparam%codename(1))
276  eq_helena(1)%codeparam%codename(1)='HELENA'
277  endif
278  if(.not.associated(eq_helena(1)%codeparam%codeversion)) then
279  allocate(eq_helena(1)%codeparam%codeversion(1))
280  eq_helena(1)%codeparam%codeversion(1)=version
281  endif
282  if(.not.associated(eq_helena(1)%codeparam%parameters)) then
283  allocate(eq_helena(1)%codeparam%parameters(1))
284  eq_helena(1)%codeparam%parameters(1)='NONE'
285  endif
286  if(.not.associated(eq_helena(1)%codeparam%output_diag)) then
287  allocate(eq_helena(1)%codeparam%output_diag(1))
288  eq_helena(1)%codeparam%output_diag(1)='NONE'
289  endif
290  eq_helena(1)%codeparam%output_flag=0
291  write(*,*) 'Calling euitm_put_non_timed for HELENA'
292  call euitm_put_non_timed(idx,"equilibrium",eq_helena(1))
293  write(*,*) 'Calling euitm_put_slice for HELENA'
294  call euitm_put_slice(idx,"equilibrium",eq_helena(1))
295  write(*,*) 'Calling euitm_closee for HELENA'
296  call euitm_close(idx)
297  write(*,*) 'Finished with euitm UAL for HELENA'
298 ! call perfoff
299 #endif
300 
301  call euitm_deallocate(eq_helena)
302 
303 #endif
304 
305  else if (eq_code.eq.2) then
306 
307 ! call helena21 and write results to 'eq_helena21.out'
308 
309 ! call perfon ('helena21')
310  call helena21itm(eq_new, eq_helena)
311 ! call perfoff
312  write(*,*) 'CODE: HELENA21'
313  call output_eq_info(eq_helena(1))
314  write(*,*)
315 
316  call open_write_file(10, 'eq_helena21.out.cpo')
317  call write_cpo(eq_helena(1), 'equilibrium')
318  call close_write_file
319  open(10,file='eq_helena21.out')
320  ! write(10,'(''# rho pressure jparallel jphi volume q rhovol psi phi vprime gm1 gm2 gm3 gm4 gm5 gm6 gm7'')')
321  write(10,'(''# pressure jphi volume q psi phi vprime gm1 gm2 gm3 gm4 gm5 gm6 gm7 rho_vol pprime ffprime F_dia area'')')
322  do irho=1,size(eq_helena(1)%profiles_1d%psi)
323  write(10,'(30f20.6)') &
324  ! EQ_HELENA(1)%profiles_1d%rho_tor(irho), &
325  eq_helena(1)%profiles_1d%pressure(irho), &
326  ! EQ_HELENA(1)%profiles_1d%jparallel(irho), &
327  eq_helena(1)%profiles_1d%jphi(irho), &
328  eq_helena(1)%profiles_1d%volume(irho), &
329  eq_helena(1)%profiles_1d%q(irho), &
330  eq_helena(1)%profiles_1d%psi(irho), &
331  eq_helena(1)%profiles_1d%phi(irho), &
332  eq_helena(1)%profiles_1d%vprime(irho), &
333  eq_helena(1)%profiles_1d%gm1(irho), &
334  eq_helena(1)%profiles_1d%gm2(irho), &
335  eq_helena(1)%profiles_1d%gm3(irho), &
336  eq_helena(1)%profiles_1d%gm4(irho), &
337  eq_helena(1)%profiles_1d%gm5(irho), &
338  eq_helena(1)%profiles_1d%gm6(irho), &
339  eq_helena(1)%profiles_1d%gm7(irho), &
340  eq_helena(1)%profiles_1d%rho_vol(irho), &
341  eq_helena(1)%profiles_1d%pprime(irho), &
342  eq_helena(1)%profiles_1d%ffprime(irho), &
343  eq_helena(1)%profiles_1d%f_dia(irho), &
344  eq_helena(1)%profiles_1d%area(irho)
345  enddo
346  close(10)
347 
348 #ifdef UAL
349 ! for HELENA21, run_no = 2
350 ! call perfon ('ualw3')
351  write(*,*) 'Calling euitm_create for HELENA21'
352  call euitm_create('euitm',shot_no,2,0,0,idx)
353  write(*,*) eq_helena(1)%time
354  eq_helena(1)%time=0.0_r8
355  if(.not.associated(eq_helena(1)%codeparam%codename)) then
356  allocate(eq_helena(1)%codeparam%codename(1))
357  eq_helena(1)%codeparam%codename(1)='HELENA21'
358  endif
359  if(.not.associated(eq_helena(1)%codeparam%codeversion)) then
360  allocate(eq_helena(1)%codeparam%codeversion(1))
361  eq_helena(1)%codeparam%codeversion(1)=version
362  endif
363  if(.not.associated(eq_helena(1)%codeparam%parameters)) then
364  allocate(eq_helena(1)%codeparam%parameters(1))
365  eq_helena(1)%codeparam%parameters(1)='NONE'
366  endif
367  if(.not.associated(eq_helena(1)%codeparam%output_diag)) then
368  allocate(eq_helena(1)%codeparam%output_diag(1))
369  eq_helena(1)%codeparam%output_diag(1)='NONE'
370  endif
371  eq_helena(1)%codeparam%output_flag=0
372  write(*,*) 'Calling euitm_put_non_timed for HELENA21'
373  call euitm_put_non_timed(idx,"equilibrium",eq_helena(1))
374  write(*,*) 'Calling euitm_put_slice for HELENA21'
375  call euitm_put_slice(idx,"equilibrium",eq_helena(1))
376  write(*,*) 'Calling euitm_close for HELENA21'
377  call euitm_close(idx)
378  write(*,*) 'Finished with euitm UAL for HELENA21'
379 ! call perfoff
380 #endif
381 
382  call euitm_deallocate(eq_helena)
383 
384  else if (eq_code.eq.3) then
385 
386 #ifdef GOT_E3EMEQ
387 
388 ! call E3EMEQ and then write results to 'eq_emeq.out'
389 
390 !dpc
391  eq_new(1)%global_param%i_plasma=5000.0_r8
392 !dpc
393 !! EQ_NEW(1)%profiles_1d%jparallel=-EQ_NEW(1)%profiles_1d%jparallel
394 ! call perfon ('emeq')
395  call emeq_e3m_wrapper(eq_new, eq_e3astr)
396 ! call perfoff
397 !! EQ_NEW(1)%profiles_1d%jparallel=-EQ_NEW(1)%profiles_1d%jparallel
398  write(*,*) 'CODE: E3EMEQ'
399  call output_eq_info(eq_e3astr(1))
400 
401  nrho1=size(eq_e3astr(1)%profiles_1d%rho_tor)
402 
403  call open_write_file(10, 'eq_emeq.out.cpo')
404  call write_cpo(eq_e3astr(1), 'equilibrium')
405  call close_write_file
406  open(10,file='eq_emeq.out')
407  write(10,'(''# rho pressure jparallel q gm1 gm2 gm3 gm4 gm5 gm6 gm7 volume vprime area aprime f_dia rho_vol'')')
408 
409  do irho=1,nrho1
410  write(10,'(30f20.6)') &
411  eq_e3astr(1)%profiles_1d%rho_tor(irho), &
412  eq_e3astr(1)%profiles_1d%pressure(irho), &
413  eq_e3astr(1)%profiles_1d%jparallel(irho), &
414  eq_e3astr(1)%profiles_1d%q(irho), &
415  eq_e3astr(1)%profiles_1d%gm1(irho), &
416  eq_e3astr(1)%profiles_1d%gm2(irho), &
417  eq_e3astr(1)%profiles_1d%gm3(irho), &
418  eq_e3astr(1)%profiles_1d%gm4(irho), &
419  eq_e3astr(1)%profiles_1d%gm5(irho), &
420  eq_e3astr(1)%profiles_1d%gm6(irho), &
421  eq_e3astr(1)%profiles_1d%gm7(irho), &
422  eq_e3astr(1)%profiles_1d%volume(irho), &
423  eq_e3astr(1)%profiles_1d%vprime(irho), &
424  eq_e3astr(1)%profiles_1d%area(irho), &
425  eq_e3astr(1)%profiles_1d%aprime(irho), &
426  eq_e3astr(1)%profiles_1d%f_dia(irho), &
427  eq_e3astr(1)%profiles_1d%rho_vol(irho)
428  enddo
429  close(10)
430 
431 #ifdef UAL
432 ! for E3ASTR, run_no=4
433 ! call perfon ('ualw4')
434  write(*,*) 'Calling euitm_create for E3ASTR'
435  call euitm_create('euitm',shot_no,4,0,0,idx)
436  write(*,*) eq_e3astr(1)%time
437  eq_e3astr(1)%time=0.0_r8
438  if(.not.associated(eq_e3astr(1)%codeparam%codename)) then
439  allocate(eq_e3astr(1)%codeparam%codename(1))
440  eq_e3astr(1)%codeparam%codename(1)='E3ASTR'
441  endif
442  if(.not.associated(eq_e3astr(1)%codeparam%codeversion)) then
443  allocate(eq_e3astr(1)%codeparam%codeversion(1))
444  eq_e3astr(1)%codeparam%codeversion(1)=version
445  endif
446  if(.not.associated(eq_e3astr(1)%codeparam%parameters)) then
447  allocate(eq_e3astr(1)%codeparam%parameters(1))
448  eq_e3astr(1)%codeparam%parameters(1)='NONE'
449  endif
450  if(.not.associated(eq_e3astr(1)%codeparam%output_diag)) then
451  allocate(eq_e3astr(1)%codeparam%output_diag(1))
452  eq_e3astr(1)%codeparam%output_diag(1)='NONE'
453  endif
454  eq_e3astr(1)%codeparam%output_flag=0
455  write(*,*) 'Calling euitm_put_non_timed for E3ASTR'
456  call euitm_put_non_timed(idx,"equilibrium",eq_e3astr(1))
457  write(*,*) 'Calling euitm_put_slice for E3ASTR'
458  call euitm_put_slice(idx,"equilibrium",eq_e3astr(1))
459  write(*,*) 'Calling euitm_close for E3ASTR'
460  call euitm_close(idx)
461  write(*,*) 'Finished with euitm UAL for E3ASTR'
462 ! call perfoff
463 #endif
464 
465  call euitm_deallocate(eq_e3astr)
466 
467 #endif
468 
469  else if(eq_code.eq.4) then
470 
471 #ifdef GOT_CHEASE
472 
473 ! call chease and write results to 'eq_chease.out'
474 
475  eq_new(1)%profiles_1d%jphi=-eq_new(1)%profiles_1d%jphi
476 ! call perfon ('chease')
477  call chease_wrapper(eq_new, eq_chease)
478 ! call perfoff
479  eq_new(1)%profiles_1d%jphi=-eq_new(1)%profiles_1d%jphi
480  write(*,*) 'CODE: CHEASE'
481  call output_eq_info(eq_chease(1))
482  write(*,*)
483 
484  call open_write_file(10, 'eq_chease.out.cpo')
485  call write_cpo(eq_chease(1), 'equilibrium')
486  call close_write_file
487  open(10,file='eq_chease.out')
488 ! write(10,'(''# rho pressure jparallel jphi volume q rhovol psi phi vprime gm1 gm2 gm3 gm4 gm5 gm6 gm7'')')
489  write(10,'(''# pressure jphi q psi volume rho_vol F_dia pprime ffprime area rho_tor phi vprime gm1 gm2 gm3 gm4 gm5 gm6 gm7'')')
490  do irho=1,size(eq_chease(1)%profiles_1d%psi)
491  write(10,'(30f20.6)') &
492  eq_chease(1)%profiles_1d%pressure(irho), &
493  eq_chease(1)%profiles_1d%jphi(irho), &
494  eq_chease(1)%profiles_1d%q(irho), &
495  eq_chease(1)%profiles_1d%psi(irho), &
496  eq_chease(1)%profiles_1d%volume(irho), &
497  eq_chease(1)%profiles_1d%rho_vol(irho), &
498  eq_chease(1)%profiles_1d%F_dia(irho), &
499  eq_chease(1)%profiles_1d%pprime(irho), &
500  eq_chease(1)%profiles_1d%ffprime(irho), &
501  eq_chease(1)%profiles_1d%area(irho), &
502  eq_chease(1)%profiles_1d%rho_tor(irho), &
503  eq_chease(1)%profiles_1d%phi(irho), &
504  eq_chease(1)%profiles_1d%vprime(irho), &
505  eq_chease(1)%profiles_1d%gm1(irho), &
506  eq_chease(1)%profiles_1d%gm2(irho), &
507  eq_chease(1)%profiles_1d%gm3(irho), &
508  eq_chease(1)%profiles_1d%gm4(irho), &
509  eq_chease(1)%profiles_1d%gm5(irho), &
510  eq_chease(1)%profiles_1d%gm6(irho), &
511  eq_chease(1)%profiles_1d%gm7(irho)
512 ! EQ_CHEASE(1)%profiles_1d%jparallel(irho), &
513  enddo
514  close(10)
515 
516 #ifdef UAL
517 ! for CHEASE, run_no=5
518 ! call perfon ('ualw5')
519  write(*,*) 'Calling euitm_create for CHEASE'
520  call euitm_create('euitm',shot_no,5,0,0,idx)
521  write(*,*) eq_chease(1)%time
522  eq_chease(1)%time=0.0_r8
523  if(.not.associated(eq_chease(1)%codeparam%codename)) then
524  allocate(eq_chease(1)%codeparam%codename(1))
525  eq_chease(1)%codeparam%codename(1)='CHEASE'
526  endif
527  if(.not.associated(eq_chease(1)%codeparam%codeversion)) then
528  allocate(eq_chease(1)%codeparam%codeversion(1))
529  eq_chease(1)%codeparam%codeversion(1)=version
530  endif
531  if(.not.associated(eq_chease(1)%codeparam%parameters)) then
532  allocate(eq_chease(1)%codeparam%parameters(1))
533  eq_chease(1)%codeparam%parameters(1)='NONE'
534  endif
535  if(.not.associated(eq_chease(1)%codeparam%output_diag)) then
536  allocate(eq_chease(1)%codeparam%output_diag(1))
537  eq_chease(1)%codeparam%output_diag(1)='NONE'
538  endif
539  eq_chease(1)%codeparam%output_flag=0
540  write(*,*) 'Calling euitm_put_non_timed for CHEASE'
541  call euitm_put_non_timed(idx,"equilibrium",eq_chease(1))
542  write(*,*) 'Calling euitm_put_slice for CHEASE'
543  call euitm_put_slice(idx,"equilibrium",eq_chease(1))
544  write(*,*) 'Calling euitm_closee for CHEASE'
545  call euitm_close(idx)
546  write(*,*) 'Finished with euitm UAL for CHEASE'
547 ! call perfoff
548 #endif
549 
550  call euitm_deallocate(eq_chease)
551 
552 #endif
553 
554  else if (eq_code.eq.5) then
555 
556 #ifdef GOT_SPIDER
557 
558 ! call SPIDER and then write results to 'spider.out'
559 
560 !dpc
561  eq_new(1)%global_param%i_plasma=5000.0_r8
562 !dpc
563 !! EQ_NEW(1)%profiles_1d%jparallel=-EQ_NEW(1)%profiles_1d%jparallel
564 ! call perfon ('emeq')
565  call spider_wrapper(eq_new, coreprof, eq_spider)
566 ! call perfoff
567 !! EQ_NEW(1)%profiles_1d%jparallel=-EQ_NEW(1)%profiles_1d%jparallel
568  write(*,*) 'CODE: SPIDER'
569  call output_eq_info(eq_spider(1))
570 
571  nrho1=size(eq_spider(1)%profiles_1d%rho_tor)
572 
573  call open_write_file(10, 'eq_emeq.out.cpo')
574  call write_cpo(eq_spider(1), 'equilibrium')
575  call close_write_file
576  open(10,file='eq_emeq.out')
577  write(10,'(''# rho pressure jparallel q gm1 gm2 gm3 gm4 gm5 gm6 gm7 volume vprime area aprime f_dia rho_vol'')')
578 
579  do irho=1,nrho1
580  write(10,'(30f20.6)') &
581  eq_spider(1)%profiles_1d%rho_tor(irho), &
582  eq_spider(1)%profiles_1d%pressure(irho), &
583  eq_spider(1)%profiles_1d%jparallel(irho), &
584  eq_spider(1)%profiles_1d%q(irho), &
585  eq_spider(1)%profiles_1d%gm1(irho), &
586  eq_spider(1)%profiles_1d%gm2(irho), &
587  eq_spider(1)%profiles_1d%gm3(irho), &
588  eq_spider(1)%profiles_1d%gm4(irho), &
589  eq_spider(1)%profiles_1d%gm5(irho), &
590  eq_spider(1)%profiles_1d%gm6(irho), &
591  eq_spider(1)%profiles_1d%gm7(irho), &
592  eq_spider(1)%profiles_1d%volume(irho), &
593  eq_spider(1)%profiles_1d%vprime(irho), &
594  eq_spider(1)%profiles_1d%area(irho), &
595  eq_spider(1)%profiles_1d%aprime(irho), &
596  eq_spider(1)%profiles_1d%f_dia(irho), &
597  eq_spider(1)%profiles_1d%rho_vol(irho)
598  enddo
599  close(10)
600 
601 #ifdef UAL
602 ! for SPIDER, run_no=6
603 ! call perfon ('ualw4')
604  write(*,*) 'Calling euitm_create for SPIDER'
605  call euitm_create('euitm',shot_no,6,0,0,idx)
606  write(*,*) eq_spider(1)%time
607  eq_spider(1)%time=0.0_r8
608  if(.not.associated(eq_spider(1)%codeparam%codename)) then
609  allocate(eq_spider(1)%codeparam%codename(1))
610  eq_spider(1)%codeparam%codename(1)='SPIDER'
611  endif
612  if(.not.associated(eq_spider(1)%codeparam%codeversion)) then
613  allocate(eq_spider(1)%codeparam%codeversion(1))
614  eq_spider(1)%codeparam%codeversion(1)=version
615  endif
616  if(.not.associated(eq_spider(1)%codeparam%parameters)) then
617  allocate(eq_spider(1)%codeparam%parameters(1))
618  eq_spider(1)%codeparam%parameters(1)='NONE'
619  endif
620  if(.not.associated(eq_spider(1)%codeparam%output_diag)) then
621  allocate(eq_spider(1)%codeparam%output_diag(1))
622  eq_spider(1)%codeparam%output_diag(1)='NONE'
623  endif
624  eq_spider(1)%codeparam%output_flag=0
625  write(*,*) 'Calling euitm_put_non_timed for SPIDER'
626  call euitm_put_non_timed(idx,"equilibrium",eq_spider(1))
627  write(*,*) 'Calling euitm_put_slice for SPIDER'
628  call euitm_put_slice(idx,"equilibrium",eq_spider(1))
629  write(*,*) 'Calling euitm_close for SPIDER'
630  call euitm_close(idx)
631  write(*,*) 'Finished with euitm UAL for SPIDER'
632 ! call perfoff
633 #endif
634 
635  call euitm_deallocate(eq_spider)
636 
637 #endif
638 
639  endif
640 
641  enddo
642 
643  call euitm_deallocate(eq_new)
644 ! call perfoff
645 ! call perfout ('whole')
646 
647 end program eq_test
648 
649 subroutine output_eq_info(eq_in)
650  use itm_constants
651  use euitm_routines
652  USE euitm_schemas
653 
654  implicit none
655 
656  TYPE (type_equilibrium) :: eq_in
657  integer iprof
658 
659  write(*,*) &
660  'psi ',associated(eq_in%profiles_1d%psi), &
661  'phi ',associated(eq_in%profiles_1d%phi), &
662  'pressure ',associated(eq_in%profiles_1d%pressure), &
663  'F_dia ',associated(eq_in%profiles_1d%F_dia), &
664  'pprime ',associated(eq_in%profiles_1d%pprime), &
665  'ffprime ',associated(eq_in%profiles_1d%ffprime), &
666  'jphi ',associated(eq_in%profiles_1d%jphi), &
667  'jparallel ',associated(eq_in%profiles_1d%jparallel), &
668  'q ',associated(eq_in%profiles_1d%q), &
669  'r_inboard ',associated(eq_in%profiles_1d%r_inboard), &
670  'r_outboard ',associated(eq_in%profiles_1d%r_outboard), &
671  'rho_vol ',associated(eq_in%profiles_1d%rho_vol), &
672  'rho_tor ',associated(eq_in%profiles_1d%rho_tor), &
673  'elongation ',associated(eq_in%profiles_1d%elongation), &
674  'tria_upper ',associated(eq_in%profiles_1d%tria_upper), &
675  'tria_lower ',associated(eq_in%profiles_1d%tria_lower), &
676  'volume ',associated(eq_in%profiles_1d%volume), &
677  'vprime ',associated(eq_in%profiles_1d%vprime), &
678  'area ',associated(eq_in%profiles_1d%area), &
679  'aprime ',associated(eq_in%profiles_1d%aprime), &
680  'gm1 ',associated(eq_in%profiles_1d%gm1), &
681  'gm2 ',associated(eq_in%profiles_1d%gm2), &
682  'gm3 ',associated(eq_in%profiles_1d%gm3), &
683  'gm4 ',associated(eq_in%profiles_1d%gm4), &
684  'gm5 ',associated(eq_in%profiles_1d%gm5), &
685  'gm6 ',associated(eq_in%profiles_1d%gm6), &
686  'gm7 ',associated(eq_in%profiles_1d%gm7), &
687  'ftrap ',associated(eq_in%profiles_1d%ftrap)
688 
689 
690  write(*,1000) ' AREA ',eq_in%global_param%area
691  write(*,1000) ' VOLUME ',eq_in%global_param%volume
692  write(*,1000) ' Raxis ',eq_in%global_param%mag_axis%position%r
693  write(*,1000) ' Zaxis ',eq_in%global_param%mag_axis%position%z
694  write(*,1000) ' Baxis ',eq_in%global_param%mag_axis%bphi
695  write(*,1000) ' Rgeo ',eq_in%eqgeometry%geom_axis%r
696  write(*,1000) ' Zgeo ',eq_in%eqgeometry%geom_axis%z
697  write(*,1000) ' a ',eq_in%eqgeometry%a_minor
698  write(*,1000) ' R0 ',eq_in%global_param%toroid_field%b0
699  write(*,1000) ' B0 ',eq_in%global_param%toroid_field%r0
700 
701  write(*,1000) ' beta_pol ',eq_in%global_param%beta_pol
702  write(*,1000) ' beta_tor ',eq_in%global_param%beta_tor
703  write(*,1000) 'beta_normal ',eq_in%global_param%beta_normal
704  write(*,1000) ' i_plasma ',eq_in%global_param%i_plasma
705  write(*,1000) ' li ',eq_in%global_param%li
706  write(*,1000) ' psi_ax ',eq_in%global_param%psi_ax
707  write(*,1000) ' psi_bound ',eq_in%global_param%psi_bound
708  write(*,1000) ' q_95 ',eq_in%global_param%q_95
709  write(*,1000) ' q_min ',eq_in%global_param%q_min
710  if(associated(eq_in%profiles_2d)) then
711  do iprof = 1, size(eq_in%profiles_2d)
712  if(associated(eq_in%profiles_2d(iprof)%grid_type)) then
713  write(*,1000) iprof, eq_in%profiles_2d(iprof)%grid_type(:)
714  endif
715  enddo
716  endif
717 
718 1000 format(a,1pg20.10)
719 
720 end subroutine output_eq_info
wrapper for HELENA
program eq_test
Framework for testing equilibrium codes.
Definition: eq_test.F90:10
wrapper for CHEASE
subroutine helena21itm(equil_in, equil_out)
Definition: helena21itm.f90:1
subroutine emeq_e3m_wrapper(EQUILIBRIUM_in, EQUILIBRIUM_out)
wrapper for BDSEQ
subroutine chease_wrapper(euitm_equilibrium_in, euitm_equilibrium_out)
subroutine helena_wrapper(euitm_equilibrium_in, euitm_equilibrium_out)
wrapper for SPIDER
program bdseq_wrapper
Definition: wrapper.F90:208
subroutine spider_wrapper(euitm_equilibrium_in, coreprof, euitm_equilibrium_out)
subroutine output_eq_info(eq_in)
Definition: eq_test.F90:649
subroutine euitm_close(idx)