PB3D  [2.45]
Ideal linear high-n MHD stability in 3-D
input_ops.f90
Go to the documentation of this file.
1 !------------------------------------------------------------------------------!
3 !------------------------------------------------------------------------------!
4 module input_ops
5 #include <PB3D_macros.h>
6  use str_utilities
7  use output_ops
8  use messages
9  use num_vars, only: dp, max_str_ln
10 
11  implicit none
12  private
14 
15 contains
22  integer function read_input_opts() result(ierr)
23  use num_vars, only: &
31  &bc_style, tol_norm, tol_slepc_loc => tol_slepc, max_it_slepc, &
41  use eq_vars, only: rho_0, r_0, pres_0, b_0, psi_0, t_0
43  &max_sec_x
45  use rich_ops, only: find_max_rich_lvl
48 
49  character(*), parameter :: rout_name = 'read_input_opts'
50 
51  ! local variables
52  character(len=max_str_ln) :: err_msg
53  integer :: pb3d_rich_lvl
54  integer, parameter :: max_size_tol_slepc = 100
55  real(dp) :: alpha
56  real(dp) :: tol_slepc(max_size_tol_slepc)
57  real(dp), parameter :: min_tol = 1.e-12_dp
58  real(dp), parameter :: max_tol = 1.e-3_dp
59 
60  ! input options
61  namelist /inputdata_pb3d/ min_par_x, max_par_x, alpha, min_r_sol, &
75  namelist /inputdata_post/ n_sol_plotted, n_theta_plot, n_zeta_plot, &
78  &pb3d_rich_lvl, max_it_zero, tol_zero, min_theta_plot, &
83  &n_vac_plot
84 
85  ! initialize ierr
86  ierr = 0
87 
88  if (rank.eq.0) then ! only master
89  call writo('Setting up user-provided input "'//trim(input_name)//&
90  &'"')
91  call lvl_ud(1)
92 
93  ! initialize input variables (optionally overwritten by user later)
94  call writo('Initialize all the inputs with default values')
95  call lvl_ud(1)
96 
97  ! common variables for all program styles
98  sol_n_procs = n_procs ! use equal amount of processes for solution as for rest of program
99  max_tot_mem = 6000_dp ! count with 6GB total
100  plot_size = [10,5] ! size of plot in inch
101  min_r_plot = 0._dp ! minimal plot normal range
102  max_r_plot = 1._dp ! maximal plot normal range
103  min_theta_plot = 1 ! starting from pi gives nicer plots
104  max_theta_plot = 3
105  select case(eq_style)
106  case (1) ! VMEC
107  n_theta_plot = 201 ! nr. poloidal points in plot
108  n_zeta_plot = 51 ! nr. toroidal points in plot
109  min_zeta_plot = 0 ! min. toroidal plot angle [pi]
110  max_zeta_plot = 2 ! max. toroidal plot angle [pi]
111  case (2) ! HELENA
112  n_theta_plot = 501 ! nr. poloidal points in plot
113  n_zeta_plot = 1 ! nr. toroidal points in plot
114  min_zeta_plot = 0 ! min. toroidal plot angle [pi]
115  max_zeta_plot = min_zeta_plot ! max. toroidal plot angle [pi]
116  end select
117  plot_grid_style = 0 ! normal plots on 3D geometry
118  max_nr_backtracks_hh = 20 ! standard nr. of backtracks
119  ex_plot_style = 1 ! GNUPlot
120  norm_disc_prec_eq = 3 ! cubic precision for normal discretization of equilibrium
121  norm_disc_prec_x = 3 ! cubic precision for normal discretization of perturbation
122  norm_disc_prec_sol = 3 ! cubic precision for normal discretization of solution
123 
124  ! select depending on program style
125  select case (prog_style)
126  case(1) ! PB3D
127  call default_input_pb3d()
128  case(2) ! POST
129  call default_input_post()
130  end select
131 
132  call lvl_ud(-1)
133 
134  ! read user input
135  ! select depending on program style
136  select case (prog_style)
137  case(1) ! PB3D
138  read(unit=input_i,nml=inputdata_pb3d,iostat=ierr,&
139  &iomsg=err_msg) ! read input data
140 
141  ! check input if successful read
142  if (ierr.eq.0) then ! input file succesfully read
143  call writo('Overwriting with user-provided file "'&
144  &//trim(input_name) // '"')
145 
146  call writo('Checking user-provided file')
147 
148  call lvl_ud(1)
149 
150  ! adapt MPI variables if needed
151  call adapt_mpi()
152 
153  ! adapt run-time variables if needed
154  ierr = adapt_run_pb3d()
155  chckerr('')
156 
157  ! adapt plotting variables if needed
158  ierr = adapt_plot()
159  chckerr('')
160 
161  ! adapt perturbation modes
162  ierr = adapt_x_modes()
163  chckerr('')
164 
165  ! adapt min_n_par_X if needed
166  call adapt_min_n_par_x
167 
168  ! adapt alpha if needed
169  ierr = adapt_alpha()
170  chckerr('')
171 
172  ! adapt solution grid
173  ierr = adapt_sol_grid(min_r_sol,max_r_sol,'sol')
174  chckerr('')
175 
176  ! adapt normalization variables if needed
177  ierr = adapt_normalization()
178  chckerr('')
179 
180  ! adapt Richardson variables if needed
181  call adapt_rich
182 
183  ! adapt input / output variables if needed
184  ierr = adapt_inoutput_pb3d()
185  chckerr('')
186 
187  ! adapt Householder variables if needed
188  call adapt_zero
189 
190  ! adapt tolerances if needed
191  call adapt_tols
192 
193  call lvl_ud(-1)
194  else ! cannot read input data
195  ierr = 1
196  call writo('Cannot read user-provided file "'&
197  &//trim(input_name)//'", error message:',&
198  &error=.true.)
199  call lvl_ud(1)
200  call writo('"'//trim(err_msg)//'"')
201  call lvl_ud(-1)
202  chckerr("")
203  end if
204  case(2) ! POST
205  read(unit=input_i,nml=inputdata_post,iostat=ierr,&
206  &iomsg=err_msg) ! read input data
207 
208  ! check input if successful read
209  if (ierr.eq.0) then ! input file succesfully read
210  call writo('Overwriting with user-provided file "'&
211  &//trim(input_name) // '"')
212 
213  call writo('Checking user-provided file')
214 
215  call lvl_ud(1)
216 
217  ! adapt run-time variables if needed
218  ierr = adapt_run_post()
219  chckerr('')
220 
221  ! adapt plotting variables if needed
222  ierr = adapt_plot()
223  chckerr('')
224 
225  ! adapt plot grid and save to solution variable
226  ierr = adapt_sol_grid(min_r_plot,max_r_plot,'plot')
227  chckerr('')
228 
229  ! adapt input / output variables if needed
230  ierr = adapt_inoutput_post()
231  chckerr('')
232 
233  call lvl_ud(-1)
234  else ! cannot read input data
235  ierr = 1
236  call writo('Cannot read user-provided file "'&
237  &//trim(input_name)//'", error message:',&
238  &error=.true.)
239  call lvl_ud(1)
240  call writo('"'//trim(err_msg)//'"')
241  call lvl_ud(-1)
242  chckerr("")
243  end if
244  end select
245 
246  ! user output
247  call writo('Close input file')
248  close(input_i)
249 
250  call lvl_ud(-1)
251  call writo('Input values set')
252  end if
253  contains
254  ! PB3D version
256  subroutine default_input_pb3d
257  !use num_vars, only: use_pol_flux_E
258 
259  ! concerning solution
260  n_r_sol = 100 ! 100 points in solution grid
261  min_r_sol = 0.1_dp ! minimum normal range
262  max_r_sol = 1.0_dp ! maximum normal range
263  tol_norm = 0.05 ! tolerance for normal range
264  ev_style = 1 ! slepc solver for EV problem
265 
266  ! concerning field line
267  alpha = 0._dp ! field line based in outboard
268  min_alpha = 0.0_dp ! minimum field-line label [pi]
269  max_alpha = 2.0_dp ! maximum field-line label [pi]
270  select case (eq_style)
271  case (1) ! VMEC
272  alpha_style = 2 ! multiple field-lines, single turns
273  n_alpha = 10
274  case (2) ! HELENA
275  alpha_style = 1 ! single field-line, multiple turns
276  end select
277 
278  ! concerning perturbation
279  min_n_par_x = req_min_n_par_x ! nonsensible value to check for user overwriting
280  min_par_x = 0.0_dp ! minimum parallel angle [pi]
281  max_par_x = 2.0_dp ! maximum parallel angle [pi]
282  prim_x = 20 ! main mode number of perturbation
283  min_sec_x = huge(1) ! nonsensible value to check for user overwriting
284  max_sec_x = huge(1) ! nonsensible value to check for user overwriting
285  n_mod_x = huge(1) ! nonsensible value to check for user overwriting
286  !use_pol_flux_F = use_pol_flux_E ! use same normal flux coordinate as the equilibrium
287  use_pol_flux_f = .true. ! use poloidail flux as toroidal flux is untested
288 
289  ! concerning normalization
290  rho_0 = 1e-7_dp ! for fusion, particle density of around 5E19 for detuerium
291  r_0 = huge(1._dp) ! nonsensible value to check for user overwriting
292  pres_0 = huge(1._dp) ! nonsensible value to check for user overwriting
293  psi_0 = huge(1._dp) ! nonsensible value to check for user overwriting
294  b_0 = huge(1._dp) ! nonsensible value to check for user overwriting
295  t_0 = huge(1._dp) ! nonsensible value to check for user overwriting
296 
297  ! concerning input / output
298  n_sol_requested = 1 ! request solution closest to guess
299  retain_all_sol = .false. ! don't retain faulty ones
300  plot_resonance = .false. ! do not plot the q-profile with nq-m = 0
301  plot_magn_grid = .false. ! do not plot the magnetic grid
302  plot_kappa = .false. ! do not plot the curvature
303  plot_b = .false. ! do not plot the magnetic field
304  plot_j = .false. ! do not plot the current
305  plot_flux_q = .false. ! do not plot the flux quantities
306 
307  ! concerning Richardson extrapolation
308  max_it_rich = 1 ! by default no Richardson extrapolation
309  tol_rich = 1.e-4_dp ! tolerance of 1E-4
310  rich_restart_lvl = 1 ! don't restart
311 
312  ! concerning finding zeros
313  max_it_zero = 100 ! maximum 100 iterations
314  tol_zero = 1.0e-10_dp ! very low tolerance for calculation of field lines
315 
316  ! runtime variables
317  use_normalization = .true. ! use normalization for the variables
318  norm_disc_style_sol = 2 ! left finite differences
319  magn_int_style = 1 ! trapezoidal rule
320  max_it_slepc = 1000 ! max. nr. of iterations for SLEPC
321  ev_bc = 1._dp ! use 1 as artificial EV for the Boundary Conditions
322  ev_guess = -0.3_dp ! sensible guess
323  tol_slepc = huge(1._dp) ! nonsensible value to check for user overwriting
324  rho_style = 1 ! constant pressure profile, equal to rho_0
325  u_style = 3 ! full expression for U, up to order 3
326  k_style = 1 ! perpendicular kinetic energy normalized
327  norm_style = 1 ! MISHKA normalization
328  bc_style = [1,2] ! left BC zeroed and right BC through minimization of energy
329  x_style = 2 ! fast style: mode numbers optimized in normal coordinate
330  solver_slepc_style = 1 ! Krylov-Schur
331  matrix_slepc_style = 1 ! sparse matrix storage
332  x_grid_style = huge(1) ! nonsensible value to check for user overwriting
333  max_njq_change = 0.49_dp ! maximum change of just under 0.5 for n q (pol. flux) or n iota (tor. flux)
334  end subroutine default_input_pb3d
335 
336  ! POST version
338  subroutine default_input_post()
339  ! concerning finding zeros
340  max_it_zero = 200 ! more iterations than PB3D
341  tol_zero = 1.0e-8_dp ! less relative error than PB3D
342 
343  ! runtime variables
344  plot_resonance = .true. ! plot the q-profile with nq-m = 0
345  plot_flux_q = .true. ! plot the flux quantities
346  plot_kappa = .true. ! plot the curvature
347  plot_magn_grid = .true. ! plot the magnetic grid
348  plot_b = .true. ! plot the magnetic field
349  plot_j = .true. ! plot the current
350  plot_sol_xi = .true. ! plot plasma perturbation of solution
351  plot_sol_q = .true. ! plot magnetic perturbation of solution
352  plot_e_rec = .true. ! plot energy reconstruction
353  plot_vac_pot = .true. ! plot vacuum potential
354  post_style = 1 ! process on extended plot grid
355  write(*,*) '!!!!!! min and max need to take into account normalization'
356  min_rvac_plot = 0.1*r_0 ! minimum R of vacuum plot
357  max_rvac_plot = 2*r_0 ! maximum R of vacuum plot
358  min_zvac_plot = -r_0 ! minimum R of vacuum plot
359  max_zvac_plot = r_0 ! maximum R of vacuum plot
360  n_vac_plot = [100,100] ! number of points in R and Z of vacuum plot
361 
362  ! variables concerning input / output
363  pert_mult_factor_post = 0._dp ! factor by which to XYZ is perturbed in POST
364  pb3d_rich_lvl = huge(1) ! don't restart
365 
366  ! variables concerning output
367  n_sol_plotted = -1 ! plot all solutions
368  end subroutine default_input_post
369 
370  ! Checks whether the variables concerning MPI are chosen correctly:
371  ! sol_n_procs cannot be greater than n_procs. If it is lower than 1,
372  ! the maximum amount of n_procs is used.
374  subroutine adapt_mpi()
375  use num_vars, only: n_procs
376 
377  if (sol_n_procs.lt.1) then
378  call writo('Using the maximum number of MPI processes '//&
379  &trim(i2str(n_procs))//' for SLEPC')
380  sol_n_procs = n_procs
381  end if
382  if (sol_n_procs.gt.n_procs) then
383  call writo('Cannot use more than '//trim(i2str(n_procs))//&
384  &' MPI processes for SLEPC',warning=.true.)
385  sol_n_procs = n_procs
386  end if
387  end subroutine adapt_mpi
388 
389  ! Checks whether the variables concerning run-time are chosen correctly:
390  ! rho_style has to be 1 (constant rho = rho_0),
391  ! matrix_SLEPC_style has to be 0 or 1,
392  ! max_it_SLEPC has to be at least 1,
393  ! magnetic integral style has to be 1..2,
394  ! perturbation grid style has to be 1..2,
395  ! tensorial perturbation variables interpolation style has to be 1..2
396  ! (for style 2, a warning should be displayed)
397  ! for HELENA (eq_style 1), only poloidal flux can be used.
398  ! norm_disc_prec_eq has to be 3 because cubics are the maximum order
399  ! for splines and derivatives of order 3 are needed for equilibrium
400  ! quantities.
401  ! The other normal discretization precisions can be lower, but not
402  ! higher.
404  integer function adapt_run_pb3d() result(ierr)
405  use num_vars, only: eq_style
406 
407  character(*), parameter :: rout_name = 'adapt_run_PB3D'
408 
409  ! initialize ierr
410  ierr = 0
411 
412  if (rho_style.ne.1) then
413  ierr = 1
414  err_msg = 'rho_style has to be 1 (constant)'
415  chckerr(err_msg)
416  end if
417  if (eq_style.eq.2 .and. .not.use_normalization) then
418  use_normalization = .true.
419  call writo('normalization is always used for HELENA',&
420  &warning=.true.)
421  end if
422  if (matrix_slepc_style.lt.1 .or. matrix_slepc_style.gt.2) then
423  ierr = 1
424  err_msg = 'matrix_SLEPC_style has to be 1 (sparse) or 2 (shell)'
425  chckerr(err_msg)
426  end if
427  if (max_it_slepc.le.0) then
428  max_it_slepc = 1
429  call writo('max_it_SLEPC has to be positive and is set to 1',&
430  &warning=.true.)
431  end if
432  if (magn_int_style.lt.1 .or. magn_int_style.gt.2) then
433  ierr = 1
434  err_msg = 'magn_int_style has to be 1 (trapezoidal) or 2 &
435  &(Simpson 3/8 rule)'
436  chckerr(err_msg)
437  end if
438  if (eq_style.eq.2 .and. .not.use_pol_flux_f) then
439  use_pol_flux_f = .true.
440  call writo('can only use poloidal flux for HELENA',&
441  &warning=.true.)
442  end if
443  if (norm_disc_prec_eq.ne.3) then
444  ierr =1
445  err_msg = 'can only use norm_disc_prec_eq = 3 currently'
446  chckerr(err_msg)
447  end if
448  if (norm_disc_prec_x.lt.1 .or. norm_disc_prec_x.gt.3) then
449  ierr =1
450  err_msg = 'norm_disc_prec_X has to be between 1 and 3'
451  chckerr(err_msg)
452  end if
453  if (norm_disc_prec_sol.lt.1 .or. norm_disc_prec_sol.gt.3) then
454  ierr = 1
455  err_msg = 'norm_disc_prec_sol has to be between 1 and 3'
456  chckerr(err_msg)
457  end if
458  end function adapt_run_pb3d
459 
460  ! Checks whether the variables concerning run-time are chosen correctly:
461  ! POST_style should be 1 (plotting on extended plot grid) or 2
462  ! (plotting on field-aligned grid also used in PB3D).
464  integer function adapt_run_post() result(ierr)
465  character(*), parameter :: rout_name = 'adapt_run_POST'
466 
467  ! initialize ierr
468  ierr = 0
469 
470  if (post_style.lt.1 .or. post_style.gt.2) then
471  ierr = 1
472  err_msg = 'POST_style has to be 1 (extended grid) or 2 &
473  &(field-aligned grid)'
474  chckerr(err_msg)
475  end if
476  end function adapt_run_post
477 
478  ! Checks whether the variables concerning output are chosen correctly:
479  ! PB3D_rich_lvl can be at most the maximally found level
480  ! POST_output_full is true if non-flux quantities are plot
481  ! POST_output_sol is true if solution quantities are plot
483  integer function adapt_inoutput_post() result(ierr)
484  use num_vars, only: compare_tor_pos
485 
486  character(*), parameter :: rout_name = 'adapt_inoutput_POST'
487 
488  ! local variables
489  integer :: max_pb3d_rich_lvl ! maximum Richardson level
490 
491  ! initialize ierr
492  ierr = 0
493 
494  ! find maximum levl
495  ierr = find_max_rich_lvl(pb3d_rich_lvl,max_pb3d_rich_lvl)
496  chckerr('')
497  if (max_pb3d_rich_lvl.le.0) then
498  call writo('No suitable Richardson level found, so only &
499  &equilibrium output will be done.',alert=.true.)
500  plot_sol_xi = .false.
501  plot_sol_q = .false.
502  plot_e_rec = .false.
503  pb3d_rich_lvl = 1
504  else
505  if (pb3d_rich_lvl.eq.huge(1)) then
506  pb3d_rich_lvl = max_pb3d_rich_lvl
507  call writo('Treating maximum Richardson level found: '//&
508  &trim(i2str(max_pb3d_rich_lvl)))
509  else if (pb3d_rich_lvl.le.max_pb3d_rich_lvl) then
510  call writo('Treating Richardson level requested: '//&
511  &trim(i2str(pb3d_rich_lvl)))
512  if (pb3d_rich_lvl.lt.max_pb3d_rich_lvl) then
513  call lvl_ud(1)
514  call writo('Less than maximum Richardson level: '//&
515  &trim(i2str(max_pb3d_rich_lvl)),warning=.true.)
516  call lvl_ud(-1)
517  end if
518  else if (pb3d_rich_lvl.gt.max_pb3d_rich_lvl) then
519  ierr = 1
520  err_msg = 'PB3D_rich_lvl cannot be higher than '//&
521  &trim(i2str(max_pb3d_rich_lvl))
522  chckerr(err_msg)
523  else if (pb3d_rich_lvl.lt.1) then
524  call writo('Richardson level requested '//&
525  &trim(i2str(pb3d_rich_lvl))//', so only equilibrium &
526  &output will be done.',alert=.true.)
527  plot_sol_xi = .false.
528  plot_sol_q = .false.
529  plot_e_rec = .false.
530  pb3d_rich_lvl = 1
531  end if
532  end if
533  max_it_rich = pb3d_rich_lvl
534  rich_lvl = pb3d_rich_lvl
535 
536  ! set POST_output_full and POST_output_sol
537  post_output_sol = plot_sol_xi .or. plot_sol_q .or. plot_e_rec
538  post_output_full = post_output_sol .or. plot_b .or. plot_j .or. &
539  &plot_kappa .or. plot_vac_pot .or. compare_tor_pos
540  end function adapt_inoutput_post
541 
542  ! Checks whether the variables concerning output are chosen correctly:
543  ! n_sol_requested has to be at least one,
544  ! rich_restart_lvl can be at most one more than the maximally found
545  ! level, nor can it be higher than max_it_rich.
547  integer function adapt_inoutput_pb3d() result(ierr)
548  character(*), parameter :: rout_name = 'adapt_inoutput_PB3D'
549 
550  ! initialize ierr
551  ierr = 0
552 
553  if (n_sol_requested.lt.1) then
554  n_sol_requested = 1
555  call writo('n_sol_requested has been increased to '&
556  &//trim(i2str(n_sol_requested)),warning=.true.)
557  end if
558  if (rich_restart_lvl.lt.1) then
559  call writo('rich_restart_lvl was '//&
560  &trim(i2str(rich_restart_lvl))//&
561  &' but cannot be lower than 1, so it was reset to 1',&
562  &warning=.true.)
563  rich_restart_lvl = 1
564  end if
565  if (rich_restart_lvl.gt.max_it_rich) then
566  ierr = 1
567  err_msg = 'rich_restart_lvl not within 1..max_it_rich = 1..'//&
568  &trim(i2str(max_it_rich))
569  chckerr(err_msg)
570  end if
571  if (rich_restart_lvl.gt.1) then
572  ierr = find_max_rich_lvl(pb3d_rich_lvl,pb3d_rich_lvl)
573  chckerr('')
574  if (rich_restart_lvl.gt.pb3d_rich_lvl+1) then
575  ierr = 1
576  if (pb3d_rich_lvl.gt.0) then
577  err_msg = 'The highest Richardson level found was '//&
578  &trim(i2str(pb3d_rich_lvl))
579  else
580  err_msg = 'No Richardson level found'
581  end if
582  chckerr(err_msg)
583  end if
584  end if
585  end function adapt_inoutput_pb3d
586 
587  ! Checks whether the variables concerning plotting are chosen correctly:
588  ! n_theta_plot and n_zeta_plot have to be positive,
589  ! if n_theta_plot or n_zeta_plot is equal to 1, the minimum value is
590  ! chosen.
591  ! if comparing different toroidal positions, only 2 values can be
592  ! used.
593  ! ex_plot_style should be 1 (GNUPlot) or 2 (Bokeh / Mayavi)
594  ! if toroidal positions are compared, theta_plot limits should contain
595  ! 0..2pi and max_r_plot should be 1.
597  integer function adapt_plot() result(ierr)
598  use num_vars, only: compare_tor_pos
599 
600  character(*), parameter :: rout_name = 'adapt_plot'
601 
602  ! local variables
603  real(dp), parameter :: tol = 1.e-6_dp ! tolerance for checks
604  character(len=max_str_ln) :: err_msg ! error message
605 
606  ! initialize ierr
607  ierr = 0
608 
609  if (n_theta_plot.lt.1) then
610  n_theta_plot = 1
611  call writo('n_theta_plot has to be positive and is set to '//&
612  &trim(i2str(n_theta_plot)),warning=.true.)
613  end if
614  if (n_zeta_plot.lt.1) then
615  n_zeta_plot = 1
616  call writo('n_zeta_plot has to be positive and is set to '//&
617  &trim(i2str(n_zeta_plot)),warning=.true.)
618  end if
619  if (n_theta_plot.eq.1 .and. &
620  &abs(max_theta_plot-min_theta_plot).gt.tol) then
621  call writo('For n_theta_plot = 1, max_theta_plot is changed to &
622  &min_theta_plot = '//trim(r2strt(min_theta_plot)),&
623  &warning=.true.)
624  max_theta_plot = min_theta_plot
625  end if
626  if (n_zeta_plot.eq.1 .and. &
627  &abs(max_zeta_plot-min_zeta_plot).gt.tol) then
628  call writo('For n_zeta_plot = 1, min_zeta_plot and &
629  &max_zeta_plot are changed to their average = '//&
630  &trim(r2strt((min_zeta_plot+max_zeta_plot)/2)),&
631  &warning=.true.)
632  max_zeta_plot = 0.5*(min_zeta_plot+max_zeta_plot)
633  min_zeta_plot = max_zeta_plot
634  end if
635  if (ex_plot_style.lt.1 .or. ex_plot_style.gt.2) then
636  call writo('external plot style should be')
637  call lvl_ud(1)
638  call writo('1: GNUPlot')
639  call writo('2: Bokeh (2-D) / Mayavi (3-D)')
640  call lvl_ud(-1)
641  call writo('Default (1: GNUPlot) chosen',warning=.true.)
642  ex_plot_style = 1
643  end if
644  if (prog_style.eq.2) then ! only for POST
645  if (compare_tor_pos) then
646  if (n_zeta_plot.ne.3) then
647  ierr = 1
648  err_msg = 'For compare_tor_pos, you can only use &
649  &n_zeta = 3 (one point in the middle)'
650  chckerr(err_msg)
651  end if
652  if (post_style.ne.1) then
653  ierr = 1
654  err_msg = 'For compare_tor_pos, you can only use &
655  &POST_style = 1 (extended grid)'
656  chckerr(err_msg)
657  end if
658  if (min_r_plot.lt.tol_zero) then
659  ierr = 1
660  err_msg = 'min_r_plot has to be greater than 0.'
661  chckerr(err_msg)
662  end if
663  if (max_r_plot.lt.1._dp) then
664  ierr = 1
665  err_msg = 'max_r_plot has to be 1.'
666  chckerr(err_msg)
667  end if
668  if (abs(min_theta_plot-0).gt.tol_zero .or. &
669  &abs(max_theta_plot-2).gt.tol_zero) then
670  ierr = 1
671  err_msg = 'theta limits of plot need to be 0..2pi'
672  chckerr(err_msg)
673  end if
674  end if
675  if (min_rvac_plot.lt.0) then
676  call writo('The minimum value for Rvac_plot has to be &
677  &greater than zero',warning=.true.)
678  call lvl_ud(1)
679  min_rvac_plot = 0.1*r_0
680  call writo('It has been set to '//&
681  &trim(r2str(min_rvac_plot)))
682  call lvl_ud(-1)
683  end if
684  if (max_rvac_plot.lt.min_rvac_plot) then
685  ierr = 1
686  err_msg = 'The maximum value for Rvac_plot has to be &
687  &greater than the minimum value'
688  chckerr(err_msg)
689  end if
690  if (eq_style.eq.1 .and. plot_vac_pot) then
691  plot_vac_pot = .false.
692  call writo('Not plotting vacuum for VMEC yet',&
693  &warning=.true.)
694  end if
695  end if
696  end function adapt_plot
697 
698  ! Checks whether variables concerning perturbation modes are correct:
699  ! the absolute value of prim_X has to be at least min_nm_X,
700  ! deduces the X style,
701  ! min_sec_X and max_sec_X, pertaining to X style 1 and n_sec_X,
702  ! pertaining to X style 2 cannot be both given,
703  ! if X style 1 (prescribed), max_sec_X has to be greater than
704  ! min_sec_X in absolute value,
705  ! if X style 2 (fast), n_mod_X has to be a number greater than 0,
706  ! for every X style, some checks are made and some variables set.
708  integer function adapt_x_modes() result(ierr)
709  use x_vars, only: min_nm_x
710 
711  character(*), parameter :: rout_name = 'adapt_X_modes'
712 
713  ! local variables
714  character(len=max_str_ln) :: err_msg ! error message
715 
716  ! initialize ierr
717  ierr = 0
718 
719  ! check prim_X: common for all X styles
720  if (abs(prim_x).lt.min_nm_x) then
721  ierr = 1
722  err_msg = 'prim_X should be at least '//trim(i2str(min_nm_x))//&
723  &', and probably a bit higher still'
724  chckerr(err_msg)
725  end if
726 
727  ! set X_style
728  if (min_sec_x.ne.huge(1) .and. min_sec_x.ne.huge(1)) then ! min_sec_X and max_sec_X prescribed
729  if (n_mod_x.eq.huge(1)) then ! n_mod_X not prescribed
730  x_style = 1 ! prescribed style
731  else ! n_mod_X also prescribed
732  ierr = 1
733  err_msg = 'Cannot prescribe both min and max of secondary &
734  &X and nr. of modes.'
735  chckerr(err_msg)
736  end if
737  else ! min_sec_X and max_sec_X not prescribed
738  if (n_mod_x.eq.huge(1)) then ! n_mod_X not prescribed either: default
739  x_style = 2 ! X style 2: fast
740  n_mod_x = 10 ! 10 resonating modes at every place
741  else ! but n_mod_X prescribed
742  x_style = 2
743  end if
744  end if
745 
746  ! set X_grid_style if not overwritten by user
747  if (x_grid_style.ge.huge(1)) then
748  select case (x_style)
749  case (1) ! prescribed
750  x_grid_style = 1 ! equilibrium
751  case (2) ! fast
752  x_grid_style = 3 ! enriched
753  end select
754  end if
755 
756  ! check X_grid style
757  if (x_grid_style.lt.1 .or. x_grid_style.gt.3) then
758  ierr = 1
759  err_msg = 'X_grid_style has to be 1 (equilibrium), 2 &
760  &(solution) or 3 (enriched)'
761  chckerr(err_msg)
762  end if
763 
764  ! checks for each X style
765  select case (x_style)
766  case (1) ! prescribed
767  ! check min_sec_X, max_sec_X
768  if (abs(max_sec_x).lt.abs(min_sec_x)) then
769  ierr = 1
770  err_msg = 'max_sec_X has to be larger or equal to &
771  &min_sec_X in absolute value'
772  chckerr(err_msg)
773  end if
774  if (max_sec_x.ge.huge(1)) then
775  ierr = 1
776  err_msg = 'max_sec_X not set'
777  chckerr(err_msg)
778  end if
779  if (min_sec_x.ge.huge(1)) then
780  ierr = 1
781  err_msg = 'min_sec_X not set'
782  chckerr(err_msg)
783  end if
784  if (x_grid_style.eq.3) then
785  ierr = 1
786  err_msg = 'X_grid_style 3 (enriched) does not make &
787  &sense for X_style 1 (prescribed)'
788  chckerr(err_msg)
789  end if
790 
791  ! set n_mod_X
792  n_mod_x = abs(max_sec_x)-(min_sec_x)+1
793  case (2) ! fast
794  if (n_mod_x.lt.1) then
795  ierr = 1
796  err_msg = 'the number of resonating modes has to be at &
797  &least one'
798  chckerr(err_msg)
799  end if
800  end select
801  end function adapt_x_modes
802 
803  ! Checks whether min_n_par_X has been chosen correctly:
804  ! It should be larger than req_min_n_par_X (arbitrary) if an integral
805  ! is to be calculated. However, there are additional requirements for
806  ! the different magnetic integral styles: for style 1 min_par_X has to
807  ! be of the form 2+k and for style 2 of the form 4+3k.
809  subroutine adapt_min_n_par_x
810  ! local variables
811  integer :: fund_n_par ! fundamental interval width
812 
813  ! check req_min_n_par_X
814  if (min_n_par_x.lt.req_min_n_par_x) then
816  call writo('min_n_par_X has been increased to '//&
817  &trim(i2str(min_n_par_x))//' to have enough parallel &
818  &resolution',warning=.true.)
819  end if
820 
821  ! check for individual magnetic integral style
822  select case (magn_int_style)
823  case (1) ! Trapezoidal rule
824  fund_n_par = 1
825  case (2) ! Simpson's 3/8 rule
826  fund_n_par = 3
827  end select
828  if (mod(min_n_par_x-1,fund_n_par).ne.0) then ! there is a remainder
829  min_n_par_x = 1 + fund_n_par * &
830  &((min_n_par_x-1)/fund_n_par + 1)
831  call writo('min_n_par_X has been increased to '//&
832  &trim(i2str(min_n_par_x))//' to be a of the form '//&
833  &trim(i2str(fund_n_par+1))//' + '//&
834  &trim(i2str(fund_n_par))//&
835  &'k for magnetic integral style '//&
836  &trim(i2str(magn_int_style)),warning=.true.)
837  end if
838  end subroutine adapt_min_n_par_x
839 
840  ! Checks whether variables related to alpha are chosen correctly :
841  ! alpha_style has to be either 0 or 1
842  ! for 1, n_alpha = 1 and min/max_alpha = alpha
843  ! for 2, max_par_X - min_par_X = 2
845  integer function adapt_alpha() result(ierr)
846  character(*), parameter :: rout_name = 'adapt_alpha'
847 
848  ! local variables
849  character(len=max_str_ln) :: err_msg ! error message
850 
851  ! initialize ierr
852  ierr = 0
853 
854  ! check alpha style
855  if (alpha_style.lt.1 .or. alpha_style.gt.2) then
856  ierr = 1
857  err_msg = 'Alpha style needs to be 1 or 2 [def]'
858  chckerr(err_msg)
859  end if
860 
861  ! adapt variables pertaining to particular alpha style
862  select case (alpha_style)
863  case (1) ! single field line, multiple turns
864  if (n_alpha.ne.1) then
865  call writo('For alpha style 1, n_alpha can only be 1',&
866  &warning=.true.)
867  n_alpha = 1
868  end if
869  min_alpha = alpha
870  max_alpha = alpha
871  case (2) ! multiple field lines, single turns
872  if (n_alpha.lt.1) then
873  ierr = 1
874  err_msg = 'For alpha style 2, n_alpha needs to be &
875  &positive'
876  chckerr(err_msg)
877  end if
878  if (abs(max_par_x - min_par_x - 2) .gt. tol_zero) then
879  min_par_x = 0._dp
880  max_par_x = 2._dp
881  call writo('min/max_par_X has been set to ['//&
882  &trim(r2strt(min_par_x))//':'//&
883  &trim(r2strt(max_par_x))//']',warning=.true.)
884  end if
885  end select
886 
887  ! for HELENA, warn about axisymmetry
888  if (eq_style.eq.2) then
889  if (abs(max_par_x - min_par_x - 2) .gt. tol_zero) then
890  call writo('HELENA equilibria are axisymmetric, so only &
891  &max - min_par_X = 2 needed',warning=.true.)
892  end if
893  if (n_alpha .gt. 1) then
894  call writo('HELENA equilibria are axisymmetric, so only &
895  &n_alpha = 1 needed',warning=.true.)
896  end if
897  end if
898  end function adapt_alpha
899 
900  ! Checks whether variables concerning the solution grid are correct:
901  ! min_r should be greater than 0
902  ! max_r should be lesser than 1
903  ! min_r should be lesser than max_r
904  ! The number of solution points should be large enough for the numerical
905  ! discretization scheme, and there should be enough points per process.
907  integer function adapt_sol_grid(min_r,max_r,var_name) result(ierr)
908  use num_vars, only: n_procs
909 
910  character(*), parameter :: rout_name = 'adapt_sol_grid'
911 
912  ! input / output
913  real(dp), intent(inout) :: min_r, max_r ! min and max of norma range
914  character(len=*), intent(in) :: var_name ! name of variable
915 
916  ! local variables
917  character(len=max_str_ln) :: err_msg ! error message
918 
919  ! initialize ierr
920  ierr = 0
921 
922  ! check min_r
923  if (min_r.lt.0._dp) then
924  min_r = 0._dp
925  call writo('min_r_'//trim(var_name)//' has been increased to '&
926  &//trim(r2strt(min_r)),warning=.true.)
927  end if
928 
929  ! check if max_r is not greater than 1
930  if (max_r.gt.1._dp) then
931  max_r = 1._dp
932  call writo('max_r_'//trim(var_name)//' has been decreased to '&
933  &//trim(r2strt(max_r)),warning=.true.)
934  end if
935 
936  ! check if min_r < max_r
937  if (min_r.ge.max_r) then
938  ierr = 1
939  err_msg = 'min_r_'//trim(var_name)//' should be lower than &
940  & max_r_'//trim(var_name)
941  chckerr(err_msg)
942  end if
943 
944  ! check n_r_sol
945  if (n_r_sol.lt.6*norm_disc_prec_sol+2) then
946  n_r_sol = 6*norm_disc_prec_sol+2
947  call writo('n_r_sol has been increased to '//&
948  &trim(i2str(n_r_sol))//' to have enough points for the &
949  &normal discretization precision',warning=.true.)
950  end if
951  if (n_r_sol.lt.n_procs*(2*norm_disc_prec_sol+1)) then
952  n_r_sol = n_procs*(2*norm_disc_prec_sol+1)
953  call writo('n_r_sol has been increased to '//&
954  &trim(i2str(n_r_sol))//' to have enough points per &
955  &process',warning=.true.)
956  end if
957 
958  ! warning if too close to axis for VMEC
959  if (eq_style.eq.1) then
960  if (min_r.le.tol_zero) &
961  &call writo('For VMEC, you should not go all the way to &
962  &the magnetic axis.',warning=.true.)
963  end if
964  end function adapt_sol_grid
965 
966  ! Checks whether the variables concerning Richardson extrapolation are
967  ! correct:
968  ! max_it_rich has to be at least 1.
970  subroutine adapt_rich
971  if (max_it_rich.lt.1) then
972  max_it_rich = 1
973  call writo('max_it_rich has been increased to 1',warning=.true.)
974  end if
975  end subroutine adapt_rich
976 
977  ! Checks whether the variables concerning finding zeros are correct:
978  ! max_it_zero has to be at least 2,
980  subroutine adapt_zero
981  if (max_it_zero.lt.1) then
982  max_it_zero = 2
983  call writo('max_it_zero has been increased to 2',warning=.true.)
984  end if
985  end subroutine adapt_zero
986 
987  ! Checks whether tolerances are correct:
988  ! tol_norm needs to be 0..1,
989  ! tol_zero needs to be min_tol..max_tol,
990  ! tol_rich needs to be min_tol..max_tol,
991  ! tol_SLEPC needs to be min_tol..max_tol.
992  ! Also sets local tol_SLEPC.
994  subroutine adapt_tols
995  ! local variables
996  integer :: id ! counter
997  integer :: max_id ! maximum of user-set id
998  real(dp), parameter :: tol_slepc_def = 1.e-5_dp ! default tol_SLEPC
999 
1000  ! check tol_norm
1001  if (tol_norm.lt.0._dp) then
1002  call writo('tol_norm has been increased to 0',warning=.true.)
1003  tol_norm = 0._dp
1004  end if
1005  if (tol_norm.gt.1._dp) then
1006  call writo('tol_norm has been decreased to 1',warning=.true.)
1007  tol_norm = 1._dp
1008  end if
1009 
1010  ! check tol_rich
1011  if (tol_rich.lt.min_tol) then
1012  call writo('tol_rich has been increased to '//&
1013  &trim(r2str(min_tol)),warning=.true.)
1014  tol_rich = min_tol
1015  end if
1016  if (tol_rich.gt.max_tol) then
1017  call writo('tol_rich has been decreased to '//&
1018  &trim(r2str(max_tol)),warning=.true.)
1019  tol_rich = max_tol
1020  end if
1021 
1022  ! check tol_zero
1023  if (tol_zero.lt.min_tol) then
1024  call writo('tol_zero has been increased to '//&
1025  &trim(r2str(min_tol)),warning=.true.)
1026  tol_zero = min_tol
1027  end if
1028  if (tol_zero.gt.max_tol) then
1029  call writo('tol_zero has been decreased to '//&
1030  &trim(r2str(max_tol)),warning=.true.)
1031  tol_zero = max_tol
1032  end if
1033 
1034  ! check and setup tol_SLEPC
1035  allocate(tol_slepc_loc(max_it_rich))
1036  ! get maximum index set by user and correct
1037  max_id = 0
1038  do id = 1,max_size_tol_slepc
1039  if (tol_slepc(id).lt.huge(1._dp)) then ! was overwritten by user
1040  max_id = id
1041  if (tol_slepc(id).lt.min_tol) then
1042  call writo('tol_SLEPC('//trim(i2str(id))//&
1043  &') has been increased to '//trim(r2str(min_tol)),&
1044  &warning=.true.)
1045  tol_slepc(id) = min_tol
1046  end if
1047  if (tol_slepc(id).gt.max_tol) then
1048  call writo('tol_SLEPC('//trim(i2str(id))//&
1049  &') has been decreased to '//trim(r2str(max_tol)),&
1050  &warning=.true.)
1051  tol_slepc(id) = max_tol
1052  end if
1053  end if
1054  end do
1055  ! set local variable
1056  if (max_id.gt.0) then ! was overwritten by user
1057  if (max_id.gt.max_it_rich) then ! too many values given
1058  call writo(trim(i2str(max_id-max_it_rich))//&
1059  &' last values discarded for tol_SLEPC',warning=.true.)
1060  tol_slepc_loc = tol_slepc(1:max_it_rich)
1061  else if (max_id.eq.max_it_rich) then ! exactly the right amount of values given
1062  tol_slepc_loc = tol_slepc(1:max_it_rich)
1063  else ! not enough values given
1064  call writo('tol_SLEPC('//trim(i2str(max_id+1))//':'//&
1065  &trim(i2str(max_it_rich))//') has been set to '//&
1066  &trim(r2str(tol_slepc_def))//' [default]',&
1067  &warning=.true.)
1068  tol_slepc_loc(1:max_id) = tol_slepc(1:max_id)
1069  tol_slepc_loc(max_id+1:max_it_rich) = tol_slepc_def
1070  end if
1071  else ! was not overwritten by user
1072  tol_slepc_loc = tol_slepc_def ! all tolerances equal to default
1073  end if
1074  end subroutine adapt_tols
1075 
1076  ! Checks whether normalization variables are chosen correctly:
1077  ! rho_0 has to be positive.
1079  integer function adapt_normalization() result(ierr)
1080  character(*), parameter :: rout_name = 'adapt_normalization'
1081 
1082  ! local variables
1083  character(len=max_str_ln) :: err_msg ! error message
1084 
1085  ! initialize ierr
1086  ierr = 0
1087 
1088  if (rho_0.le.0) then
1089  ierr = 1
1090  err_msg = 'rho_0 has to be positive'
1091  chckerr(err_msg)
1092  end if
1093  end function adapt_normalization
1094  end function read_input_opts
1095 
1103  integer function read_input_eq() result(ierr)
1105  use vmec_ops, only: read_vmec
1106  use helena_ops, only: read_hel
1107  use grid_vars, only: n_r_in
1108 
1109  character(*), parameter :: rout_name = 'read_input_eq'
1110 
1111 
1112  ! initialize ierr
1113  ierr = 0
1114 
1115  ! choose which equilibrium style is being used:
1116  ! 1: VMEC
1117  ! 2: HELENA
1118  select case (eq_style)
1119  case (1) ! VMEC
1121  chckerr('')
1122  case (2) ! HELENA
1124  chckerr('')
1125  end select
1126 
1127  ! close equilibrium file
1128  close(eq_i)
1129  end function read_input_eq
1130 
1173  integer function print_output_in(data_name,remove_previous_arrs) &
1174  &result(ierr)
1182  use eq_vars, only: r_0, pres_0, b_0, psi_0, rho_0, t_0, vac_perm, &
1184  use grid_vars, only: n_r_in, n_r_eq, n_r_sol, n_alpha, min_alpha, &
1185  &max_alpha
1186  use grid_ops, only: calc_norm_range
1188  &n_mod_x
1189  use hdf5_ops, only: print_hdf5_arrs
1190  use hdf5_vars, only: dealloc_var_1d, var_1d_type, &
1192  use helena_vars, only: chi_h, flux_p_h, flux_t_h, r_h, z_h, nchi, ias, &
1194  use vmec_vars, only: is_freeb_v, mnmax_v, mpol_v, ntor_v, is_asym_v, &
1195  &gam_v, r_v_c, r_v_s, z_v_c, z_v_s, l_v_c, l_v_s, jac_v_c, &
1196  &jac_v_s, mnmax_v, mn_v, rot_t_v, q_saf_v, pres_v, flux_t_v, &
1197  &flux_p_v, nfp_v, b_v_sub_c, b_v_sub_s
1198  use helena_vars, only: h_h_11, h_h_12, h_h_33
1199 #if ldebug
1200  use vmec_vars, only: b_v_c, b_v_s, j_v_sup_int
1201 #endif
1202 
1203  character(*), parameter :: rout_name = 'print_output_in'
1204 
1205  ! input / output
1206  character(len=*), intent(in) :: data_name
1207  logical, intent(in), optional :: remove_previous_arrs
1208 
1209  ! local variables
1210  type(var_1d_type), allocatable, target :: in_1d(:) ! 1D equivalent of input variables
1211  type(var_1d_type), pointer :: in_1d_loc => null() ! local element in in_1D
1212  integer :: id ! counter
1213  integer :: in_limits(2) ! min. and max. index of input variable grid of this process
1214  integer :: loc_size ! local size
1215 
1216  ! initialize ierr
1217  ierr = 0
1218 
1219  ! user output
1220  call writo('Write input variables to output file')
1221  call lvl_ud(1)
1222 
1223  ! plot modes if requested and VMEC
1224  if (plot_vmec_modes) then
1225  call writo('Plotting decay of VMEC modes')
1226  call print_ex_2d('R_V_c','R_V_c',&
1227  &log10(maxval(abs(r_v_c(:,:,0)),2)),draw=.false.)
1228  call draw_ex(['R_V_c'],'R_V_c',1,1,.false.,ex_plot_style=1)
1229  call print_ex_2d('R_V_s','R_V_s',&
1230  &log10(maxval(abs(r_v_s(:,:,0)),2)),draw=.false.)
1231  call draw_ex(['R_V_s'],'R_V_s',1,1,.false.,ex_plot_style=1)
1232  call print_ex_2d('Z_V_c','Z_V_c',&
1233  &log10(maxval(abs(z_v_c(:,:,0)),2)),draw=.false.)
1234  call draw_ex(['Z_V_c'],'Z_V_c',1,1,.false.,ex_plot_style=1)
1235  call print_ex_2d('Z_V_s','Z_V_s',&
1236  &log10(maxval(abs(z_v_s(:,:,0)),2)),draw=.false.)
1237  call draw_ex(['Z_V_s'],'Z_V_s',1,1,.false.,ex_plot_style=1)
1238  call print_ex_2d('L_V_c','L_V_c',&
1239  &log10(maxval(abs(l_v_c(:,:,0)),2)),draw=.false.)
1240  call draw_ex(['L_V_c'],'L_V_c',1,1,.false.,ex_plot_style=1)
1241  call print_ex_2d('L_V_s','L_V_s',&
1242  &log10(maxval(abs(l_v_s(:,:,0)),2)),draw=.false.)
1243  call draw_ex(['L_V_s'],'L_V_s',1,1,.false.,ex_plot_style=1)
1244  call print_ex_2d('jac_V_c','jac_V_c',&
1245  &log10(maxval(abs(jac_v_c(:,:,0)),2)),draw=.false.)
1246  call draw_ex(['jac_V_c'],'jac_V_c',1,1,.false.,ex_plot_style=1)
1247  call print_ex_2d('jac_V_s','jac_V_s',&
1248  &log10(maxval(abs(jac_v_s(:,:,0)),2)),draw=.false.)
1249  call draw_ex(['jac_V_s'],'jac_V_s',1,1,.false.,ex_plot_style=1)
1250 #if ldebug
1251  call print_ex_2d('B_V_c','B_V_c',&
1252  &log10(maxval(abs(b_v_c),2)),draw=.false.)
1253  call draw_ex(['B_V_c'],'B_V_c',1,1,.false.,ex_plot_style=1)
1254  call print_ex_2d('B_V_s','B_V_s',&
1255  &log10(maxval(abs(b_v_s),2)),draw=.false.)
1256  call draw_ex(['B_V_s'],'B_V_s',1,1,.false.,ex_plot_style=1)
1257 #endif
1258  end if
1259 
1260  ! calculate limits of input range
1261  ierr = calc_norm_range('PB3D_in',in_limits=in_limits)
1262  chckerr('')
1263  n_r_eq = in_limits(2)-in_limits(1)+1
1264  if (in_limits(1).gt.1 .or. in_limits(2).lt.n_r_in) then
1265  call writo('Only the normal range '//trim(i2str(in_limits(1)))//&
1266  &'..'//trim(i2str(in_limits(2)))//' is written')
1267  else
1268  call writo('The entire normal range '//trim(i2str(in_limits(1)))//&
1269  &'..'//trim(i2str(in_limits(2)))//' is written')
1270  end if
1271  call writo('(relevant to solution range '//&
1272  &trim(r2strt(min_r_sol))//'..'//trim(r2strt(max_r_sol))//')')
1273 
1274  ! set up 1D equivalents of input variables
1275  allocate(in_1d(max_dim_var_1d))
1276 
1277  ! Set up common variables in_1D
1278  id = 1
1279 
1280  ! misc_in
1281  in_1d_loc => in_1d(id); id = id+1
1282  in_1d_loc%var_name = 'misc_in'
1283  allocate(in_1d_loc%tot_i_min(1),in_1d_loc%tot_i_max(1))
1284  allocate(in_1d_loc%loc_i_min(1),in_1d_loc%loc_i_max(1))
1285  in_1d_loc%loc_i_min = [1]
1286  in_1d_loc%loc_i_max = [20]
1287  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1288  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1289  allocate(in_1d_loc%p(20))
1290  in_1d_loc%p = [prog_version,eq_style*1._dp,rho_style*1._dp,&
1291  &r_0,pres_0,b_0,psi_0,rho_0,t_0,vac_perm,-1._dp,-1._dp,&
1292  &-1._dp,norm_disc_prec_eq*1._dp,n_r_in*1._dp,n_r_eq*1._dp,&
1293  &n_r_sol*1._dp,max_flux_e,max_flux_f,-1._dp]
1294  if (use_pol_flux_e) in_1d_loc%p(11) = 1._dp
1295  if (use_pol_flux_f) in_1d_loc%p(12) = 1._dp
1296  if (use_normalization) in_1d_loc%p(13) = 1._dp
1297  if (debug_version) in_1d_loc%p(20) = 1._dp
1298 
1299  ! misc_in_V or misc_in_H, depending on equilibrium style
1300  select case (eq_style)
1301  case (1) ! VMEC
1302  ! misc_in_V
1303  in_1d_loc => in_1d(id); id = id+1
1304  in_1d_loc%var_name = 'misc_in_V'
1305  allocate(in_1d_loc%tot_i_min(1),&
1306  &in_1d_loc%tot_i_max(1))
1307  allocate(in_1d_loc%loc_i_min(1),&
1308  &in_1d_loc%loc_i_max(1))
1309  in_1d_loc%loc_i_min = [1]
1310  in_1d_loc%loc_i_max = [7]
1311  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1312  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1313  allocate(in_1d_loc%p(7))
1314  in_1d_loc%p = [-1._dp,-1._dp,mnmax_v*1._dp,mpol_v*1._dp,&
1315  &ntor_v*1._dp,nfp_v*1._dp,gam_v]
1316  if (is_asym_v) in_1d_loc%p(1) = 1._dp
1317  if (is_freeb_v) in_1d_loc%p(2) = 1._dp
1318 
1319  loc_size = n_r_eq*size(flux_t_v,2)
1320 
1321  ! flux_t_V
1322  in_1d_loc => in_1d(id); id = id+1
1323  in_1d_loc%var_name = 'flux_t_V'
1324  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1325  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1326  in_1d_loc%loc_i_min = [1,0]
1327  in_1d_loc%loc_i_max = [n_r_eq,size(flux_t_v,2)-1]
1328  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1329  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1330  allocate(in_1d_loc%p(loc_size))
1331  in_1d_loc%p = reshape(flux_t_v(in_limits(1):in_limits(2),:),&
1332  &[loc_size])
1333 
1334  ! flux_p_V
1335  in_1d_loc => in_1d(id); id = id+1
1336  in_1d_loc%var_name = 'flux_p_V'
1337  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1338  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1339  in_1d_loc%loc_i_min = [1,0]
1340  in_1d_loc%loc_i_max = [n_r_eq,size(flux_p_v,2)-1]
1341  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1342  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1343  allocate(in_1d_loc%p(loc_size))
1344  in_1d_loc%p = reshape(flux_p_v(in_limits(1):in_limits(2),:),&
1345  &[loc_size])
1346 
1347  loc_size = n_r_eq*size(pres_v,2)
1348 
1349  ! pres_V
1350  in_1d_loc => in_1d(id); id = id+1
1351  in_1d_loc%var_name = 'pres_V'
1352  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1353  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1354  in_1d_loc%loc_i_min = [1,0]
1355  in_1d_loc%loc_i_max = [n_r_eq,size(pres_v,2)-1]
1356  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1357  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1358  allocate(in_1d_loc%p(loc_size))
1359  in_1d_loc%p = reshape(pres_v(in_limits(1):in_limits(2),:),&
1360  &[loc_size])
1361 
1362  ! rot_t_V
1363  in_1d_loc => in_1d(id); id = id+1
1364  in_1d_loc%var_name = 'rot_t_V'
1365  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1366  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1367  in_1d_loc%loc_i_min = [1,0]
1368  in_1d_loc%loc_i_max = [n_r_eq,size(rot_t_v,2)-1]
1369  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1370  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1371  allocate(in_1d_loc%p(loc_size))
1372  in_1d_loc%p = reshape(rot_t_v(in_limits(1):in_limits(2),:),&
1373  &[loc_size])
1374 
1375  ! q_saf_V
1376  in_1d_loc => in_1d(id); id = id+1
1377  in_1d_loc%var_name = 'q_saf_V'
1378  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1379  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1380  in_1d_loc%loc_i_min = [1,0]
1381  in_1d_loc%loc_i_max = [n_r_eq,size(q_saf_v,2)-1]
1382  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1383  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1384  allocate(in_1d_loc%p(loc_size))
1385  in_1d_loc%p = reshape(q_saf_v(in_limits(1):in_limits(2),:),&
1386  &[loc_size])
1387 
1388  ! mn_V
1389  in_1d_loc => in_1d(id); id = id+1
1390  in_1d_loc%var_name = 'mn_V'
1391  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1392  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1393  in_1d_loc%loc_i_min = [1,1]
1394  in_1d_loc%loc_i_max = [mnmax_v,2]
1395  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1396  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1397  allocate(in_1d_loc%p(size(mn_v)))
1398  in_1d_loc%p = reshape(mn_v,[size(mn_v)])
1399 
1400  ! RZL_V
1401  in_1d_loc => in_1d(id); id = id+1
1402  in_1d_loc%var_name = 'RZL_V'
1403  allocate(in_1d_loc%tot_i_min(4),in_1d_loc%tot_i_max(4))
1404  allocate(in_1d_loc%loc_i_min(4),in_1d_loc%loc_i_max(4))
1405  in_1d_loc%loc_i_min = [1,1,0,1]
1406  in_1d_loc%loc_i_max = [mnmax_v,n_r_eq,size(r_v_c,3)-1,6]
1407  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1408  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1409  allocate(in_1d_loc%p(6*mnmax_v*n_r_eq*size(r_v_c,3)))
1410  in_1d_loc%p = reshape([r_v_c(:,in_limits(1):in_limits(2),:),&
1411  &r_v_s(:,in_limits(1):in_limits(2),:),&
1412  &z_v_c(:,in_limits(1):in_limits(2),:),&
1413  &z_v_s(:,in_limits(1):in_limits(2),:),&
1414  &l_v_c(:,in_limits(1):in_limits(2),:),&
1415  &l_v_s(:,in_limits(1):in_limits(2),:)],&
1416  &[6*mnmax_v*n_r_eq*size(r_v_c,3)])
1417 
1418  ! jac_V
1419  in_1d_loc => in_1d(id); id = id+1
1420  in_1d_loc%var_name = 'jac_V'
1421  allocate(in_1d_loc%tot_i_min(4),in_1d_loc%tot_i_max(4))
1422  allocate(in_1d_loc%loc_i_min(4),in_1d_loc%loc_i_max(4))
1423  in_1d_loc%loc_i_min = [1,1,0,1]
1424  in_1d_loc%loc_i_max = [mnmax_v,n_r_eq,size(jac_v_c,3)-1,2]
1425  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1426  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1427  allocate(in_1d_loc%p(2*mnmax_v*n_r_eq*size(jac_v_c,3)))
1428  in_1d_loc%p = reshape([jac_v_c(:,in_limits(1):in_limits(2),:),&
1429  &jac_v_s(:,in_limits(1):in_limits(2),:)],&
1430  &[2*mnmax_v*n_r_eq*size(jac_v_c,3)])
1431 
1432  ! B_V_sub
1433  in_1d_loc => in_1d(id); id = id+1
1434  in_1d_loc%var_name = 'B_V_sub'
1435  allocate(in_1d_loc%tot_i_min(4),in_1d_loc%tot_i_max(4))
1436  allocate(in_1d_loc%loc_i_min(4),in_1d_loc%loc_i_max(4))
1437  in_1d_loc%loc_i_min = [1,1,1,1]
1438  in_1d_loc%loc_i_max = [mnmax_v,n_r_eq,3,2]
1439  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1440  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1441  allocate(in_1d_loc%p(2*mnmax_v*n_r_eq*3))
1442  in_1d_loc%p = reshape([&
1443  &b_v_sub_c(:,in_limits(1):in_limits(2),:),&
1444  &b_v_sub_s(:,in_limits(1):in_limits(2),:)],&
1445  &[2*mnmax_v*n_r_eq*3])
1446 
1447 #if ldebug
1448  ! B_V
1449  in_1d_loc => in_1d(id); id = id+1
1450  in_1d_loc%var_name = 'B_V'
1451  allocate(in_1d_loc%tot_i_min(3),in_1d_loc%tot_i_max(3))
1452  allocate(in_1d_loc%loc_i_min(3),in_1d_loc%loc_i_max(3))
1453  in_1d_loc%loc_i_min = [1,1,1]
1454  in_1d_loc%loc_i_max = [mnmax_v,n_r_eq,2]
1455  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1456  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1457  allocate(in_1d_loc%p(2*mnmax_v*n_r_eq))
1458  in_1d_loc%p = reshape([b_v_c(:,in_limits(1):in_limits(2)),&
1459  &b_v_s(:,in_limits(1):in_limits(2))],[2*mnmax_v*n_r_eq])
1460 
1461  ! J_V_sup_int
1462  in_1d_loc => in_1d(id); id = id+1
1463  in_1d_loc%var_name = 'J_V_sup_int'
1464  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1465  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1466  in_1d_loc%loc_i_min = [1,1]
1467  in_1d_loc%loc_i_max = [n_r_eq,2]
1468  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1469  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1470  allocate(in_1d_loc%p(2*n_r_eq))
1471  in_1d_loc%p = reshape(j_v_sup_int(in_limits(1):in_limits(2),:),&
1472  &[2*n_r_eq])
1473 
1474 #endif
1475  case (2) ! HELENA
1476  ! misc_in_H
1477  in_1d_loc => in_1d(id); id = id+1
1478  in_1d_loc%var_name = 'misc_in_H'
1479  allocate(in_1d_loc%tot_i_min(1),&
1480  &in_1d_loc%tot_i_max(1))
1481  allocate(in_1d_loc%loc_i_min(1),&
1482  &in_1d_loc%loc_i_max(1))
1483  in_1d_loc%loc_i_min = [1]
1484  in_1d_loc%loc_i_max = [2]
1485  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1486  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1487  allocate(in_1d_loc%p(2))
1488  in_1d_loc%p = [ias*1._dp,nchi*1._dp]
1489 
1490  ! RZ_H
1491  in_1d_loc => in_1d(id); id = id+1
1492  in_1d_loc%var_name = 'RZ_H'
1493  allocate(in_1d_loc%tot_i_min(3),in_1d_loc%tot_i_max(3))
1494  allocate(in_1d_loc%loc_i_min(3),in_1d_loc%loc_i_max(3))
1495  in_1d_loc%loc_i_min = [1,1,1]
1496  in_1d_loc%loc_i_max = [nchi,n_r_eq,2]
1497  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1498  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1499  allocate(in_1d_loc%p(2*nchi*n_r_eq))
1500  in_1d_loc%p = reshape([r_h(:,in_limits(1):in_limits(2)),&
1501  &z_h(:,in_limits(1):in_limits(2))],[2*nchi*n_r_eq])
1502 
1503  ! chi_H
1504  in_1d_loc => in_1d(id); id = id+1
1505  in_1d_loc%var_name = 'chi_H'
1506  allocate(in_1d_loc%tot_i_min(1),in_1d_loc%tot_i_max(1))
1507  allocate(in_1d_loc%loc_i_min(1),in_1d_loc%loc_i_max(1))
1508  in_1d_loc%loc_i_min = [1]
1509  in_1d_loc%loc_i_max = [nchi]
1510  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1511  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1512  allocate(in_1d_loc%p(nchi))
1513  in_1d_loc%p = chi_h
1514 
1515  loc_size = n_r_eq*size(flux_p_h,2)
1516 
1517  ! flux_p_H
1518  in_1d_loc => in_1d(id); id = id+1
1519  in_1d_loc%var_name = 'flux_p_H'
1520  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1521  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1522  in_1d_loc%loc_i_min = [1,0]
1523  in_1d_loc%loc_i_max = [n_r_eq,size(flux_p_h,2)-1]
1524  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1525  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1526  allocate(in_1d_loc%p(loc_size))
1527  in_1d_loc%p = reshape(flux_p_h(in_limits(1):in_limits(2),:),&
1528  &[loc_size])
1529 
1530  ! flux_t_H
1531  in_1d_loc => in_1d(id); id = id+1
1532  in_1d_loc%var_name = 'flux_t_H'
1533  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1534  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1535  in_1d_loc%loc_i_min = [1,0]
1536  in_1d_loc%loc_i_max = [n_r_eq,size(flux_t_h,2)-1]
1537  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1538  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1539  allocate(in_1d_loc%p(loc_size))
1540  in_1d_loc%p = reshape(flux_t_h(in_limits(1):in_limits(2),:),&
1541  &[loc_size])
1542 
1543  ! q_saf_H
1544  in_1d_loc => in_1d(id); id = id+1
1545  in_1d_loc%var_name = 'q_saf_H'
1546  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1547  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1548  in_1d_loc%loc_i_min = [1,0]
1549  in_1d_loc%loc_i_max = [n_r_eq,size(q_saf_h,2)-1]
1550  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1551  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1552  allocate(in_1d_loc%p(loc_size))
1553  in_1d_loc%p = reshape(q_saf_h(in_limits(1):in_limits(2),:),&
1554  &[loc_size])
1555 
1556  ! rot_t_H
1557  in_1d_loc => in_1d(id); id = id+1
1558  in_1d_loc%var_name = 'rot_t_H'
1559  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1560  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1561  in_1d_loc%loc_i_min = [1,0]
1562  in_1d_loc%loc_i_max = [n_r_eq,size(rot_t_h,2)-1]
1563  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1564  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1565  allocate(in_1d_loc%p(loc_size))
1566  in_1d_loc%p = reshape(rot_t_h(in_limits(1):in_limits(2),:),&
1567  &[loc_size])
1568 
1569  ! pres_H
1570  in_1d_loc => in_1d(id); id = id+1
1571  in_1d_loc%var_name = 'pres_H'
1572  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1573  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1574  in_1d_loc%loc_i_min = [1,0]
1575  in_1d_loc%loc_i_max = [n_r_eq,size(pres_h,2)-1]
1576  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1577  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1578  allocate(in_1d_loc%p(loc_size))
1579  in_1d_loc%p = reshape(pres_h(in_limits(1):in_limits(2),:),&
1580  &[loc_size])
1581 
1582  ! RBphi_H
1583  in_1d_loc => in_1d(id); id = id+1
1584  in_1d_loc%var_name = 'RBphi_H'
1585  allocate(in_1d_loc%tot_i_min(2),in_1d_loc%tot_i_max(2))
1586  allocate(in_1d_loc%loc_i_min(2),in_1d_loc%loc_i_max(2))
1587  in_1d_loc%loc_i_min = [1,0]
1588  in_1d_loc%loc_i_max = [n_r_eq,size(rbphi_h,2)-1]
1589  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1590  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1591  allocate(in_1d_loc%p(loc_size))
1592  in_1d_loc%p = reshape(rbphi_h(in_limits(1):in_limits(2),:),&
1593  &[loc_size])
1594 
1595  ! h_H
1596  in_1d_loc => in_1d(id); id = id+1
1597  in_1d_loc%var_name = 'h_H'
1598  allocate(in_1d_loc%tot_i_min(3),in_1d_loc%tot_i_max(3))
1599  allocate(in_1d_loc%loc_i_min(3),in_1d_loc%loc_i_max(3))
1600  in_1d_loc%loc_i_min = [1,1,1]
1601  in_1d_loc%loc_i_max = [nchi,n_r_eq,3]
1602  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1603  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1604  allocate(in_1d_loc%p(3*nchi*n_r_eq))
1605  in_1d_loc%p = reshape([h_h_11(:,in_limits(1):in_limits(2)),&
1606  &h_h_12(:,in_limits(1):in_limits(2)),&
1607  &h_h_33(:,in_limits(1):in_limits(2))],&
1608  &[3*nchi*n_r_eq])
1609  end select
1610 
1611  ! misc_X
1612  in_1d_loc => in_1d(id); id = id+1
1613  in_1d_loc%var_name = 'misc_X'
1614  allocate(in_1d_loc%tot_i_min(1),in_1d_loc%tot_i_max(1))
1615  allocate(in_1d_loc%loc_i_min(1),in_1d_loc%loc_i_max(1))
1616  in_1d_loc%loc_i_min = [1]
1617  in_1d_loc%loc_i_max = [17]
1618  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1619  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1620  allocate(in_1d_loc%p(17))
1621  in_1d_loc%p = [prim_x*1._dp,n_mod_x*1._dp,min_sec_x*1._dp,&
1622  &max_sec_x*1._dp,norm_disc_prec_x*1._dp,norm_style*1._dp,&
1623  &u_style*1._dp,x_style*1._dp,matrix_slepc_style*1._dp,&
1624  &magn_int_style*1._dp,k_style*1._dp,alpha_style*1._dp,&
1625  &x_grid_style*1._dp,min_alpha,max_alpha,n_alpha*1._dp,&
1626  &max_njq_change]
1627 
1628  ! misc_sol
1629  in_1d_loc => in_1d(id); id = id+1
1630  in_1d_loc%var_name = 'misc_sol'
1631  allocate(in_1d_loc%tot_i_min(1),in_1d_loc%tot_i_max(1))
1632  allocate(in_1d_loc%loc_i_min(1),in_1d_loc%loc_i_max(1))
1633  in_1d_loc%loc_i_min = [1]
1634  in_1d_loc%loc_i_max = [8]
1635  in_1d_loc%tot_i_min = in_1d_loc%loc_i_min
1636  in_1d_loc%tot_i_max = in_1d_loc%loc_i_max
1637  allocate(in_1d_loc%p(8))
1638  in_1d_loc%p = [min_r_sol,max_r_sol,norm_disc_prec_sol*1._dp,&
1639  &norm_disc_style_sol*1._dp,bc_style(1)*1._dp,bc_style(2)*1._dp,&
1640  &ev_style*1._dp,ev_bc]
1641 
1642  ! write
1643  ierr = print_hdf5_arrs(in_1d(1:id-1),pb3d_name,trim(data_name),&
1644  &remove_previous_arrs=remove_previous_arrs)
1645  chckerr('')
1646 
1647  ! clean up
1648  call dealloc_var_1d(in_1d)
1649  nullify(in_1d_loc)
1650 
1651  ! user output
1652  call lvl_ud(-1)
1653  end function print_output_in
1654 end module input_ops
num_vars::plot_e_rec
logical, public plot_e_rec
whether to plot energy reconstruction in POST
Definition: num_vars.f90:111
num_vars::x_grid_style
integer, public x_grid_style
style for normal component of X grid (1: eq, 2: sol, 3: enriched)
Definition: num_vars.f90:99
num_vars::ex_plot_style
integer, public ex_plot_style
external plot style (1: GNUPlot, 2: Bokeh for 2D, Mayavi for 3D)
Definition: num_vars.f90:175
num_vars::u_style
integer, public u_style
style for calculation of U (1: ord.2, 2: ord.1, 1: ord.0)
Definition: num_vars.f90:91
num_vars::plot_sol_xi
logical, public plot_sol_xi
whether to plot plasma perturbation of solution in POST
Definition: num_vars.f90:108
num_vars::plot_resonance
logical, public plot_resonance
whether to plot the q-profile or iota-profile with resonances
Definition: num_vars.f90:102
hdf5_vars::dealloc_var_1d
Deallocates 1D variable.
Definition: HDF5_vars.f90:68
vmec_ops
Operations that concern the output of VMEC.
Definition: VMEC_ops.f90:4
helena_vars::rbphi_h
real(dp), dimension(:,:), allocatable, public rbphi_h
Definition: HELENA_vars.f90:29
num_vars::dp
integer, parameter, public dp
double precision
Definition: num_vars.f90:46
eq_vars
Variables that have to do with equilibrium quantities and the grid used in the calculations:
Definition: eq_vars.f90:27
num_vars::input_i
integer, parameter, public input_i
file number of input file
Definition: num_vars.f90:182
num_vars::n_sol_plotted
integer, dimension(4), public n_sol_plotted
how many solutions to be plot (first unstable, last unstable, first stable, last stable)
Definition: num_vars.f90:171
num_vars::use_normalization
logical, public use_normalization
whether to use normalization or not
Definition: num_vars.f90:115
eq_vars::vac_perm
real(dp), public vac_perm
either usual mu_0 (default) or normalized
Definition: eq_vars.f90:48
num_vars::bc_style
integer, dimension(2), public bc_style
style for BC left and right
Definition: num_vars.f90:94
num_vars
Numerical variables used by most other modules.
Definition: num_vars.f90:4
num_vars::input_name
character(len=max_str_ln), public input_name
will hold the full name of the input file
Definition: num_vars.f90:174
vmec_vars::mn_v
integer, dimension(:,:), allocatable, public mn_v
m and n of modes
Definition: VMEC_vars.f90:32
num_vars::max_str_ln
integer, parameter, public max_str_ln
maximum length of strings
Definition: num_vars.f90:50
num_vars::max_it_zero
integer, public max_it_zero
maximum number of iterations to find zeros
Definition: num_vars.f90:131
helena_vars::rot_t_h
real(dp), dimension(:,:), allocatable, public rot_t_h
rotational transform
Definition: HELENA_vars.f90:28
num_vars::plot_sol_q
logical, public plot_sol_q
whether to plot magnetic perturbation of solution in POST
Definition: num_vars.f90:109
rich_vars
Variables concerning Richardson extrapolation.
Definition: rich_vars.f90:4
num_vars::norm_disc_prec_eq
integer, public norm_disc_prec_eq
precision for normal discretization for equilibrium
Definition: num_vars.f90:120
helena_vars::q_saf_h
real(dp), dimension(:,:), allocatable, public q_saf_h
safety factor
Definition: HELENA_vars.f90:27
str_utilities::i2str
elemental character(len=max_str_ln) function, public i2str(k)
Convert an integer to string.
Definition: str_utilities.f90:18
num_vars::retain_all_sol
logical, public retain_all_sol
retain also faulty solutions
Definition: num_vars.f90:152
hdf5_ops
Operations on HDF5 and XDMF variables.
Definition: HDF5_ops.f90:27
grid_vars::max_par_x
real(dp), public max_par_x
max. of parallel coordinate [ ] in field-aligned grid
Definition: grid_vars.f90:25
num_vars::plot_size
integer, dimension(2), public plot_size
size of plot in inches
Definition: num_vars.f90:172
num_vars::n_sol_requested
integer, public n_sol_requested
how many solutions requested
Definition: num_vars.f90:170
helena_vars::flux_p_h
real(dp), dimension(:,:), allocatable, public flux_p_h
poloidal flux
Definition: HELENA_vars.f90:24
num_vars::max_njq_change
real(dp), public max_njq_change
maximum change of prim. mode number times saf. fac. / rot. transf. when using X_style 2 (fast)
Definition: num_vars.f90:119
num_vars::n_zeta_plot
integer, public n_zeta_plot
nr. of toroidal points in plot
Definition: num_vars.f90:157
vmec_vars::b_v_c
real(dp), dimension(:,:), allocatable, public b_v_c
Coeff. of magnitude of B in sine series (HM and FM)
Definition: VMEC_vars.f90:50
x_vars::min_r_sol
real(dp), public min_r_sol
min. normal range for pert.
Definition: X_vars.f90:135
vmec_vars::j_v_sup_int
real(dp), dimension(:,:), allocatable, public j_v_sup_int
Integrated poloidal and toroidal current (FM)
Definition: VMEC_vars.f90:52
num_vars::post_output_sol
logical, public post_output_sol
POST has outputs of solution.
Definition: num_vars.f90:147
num_vars::min_zvac_plot
real(dp), public min_zvac_plot
min. of Z in which to plot vacuum
Definition: num_vars.f90:167
eq_vars::t_0
real(dp), public t_0
derived normalization constant for nondimensionalization
Definition: eq_vars.f90:47
num_vars::max_theta_plot
real(dp), public max_theta_plot
max. of theta_plot
Definition: num_vars.f90:160
num_vars::n_procs
integer, public n_procs
nr. of MPI processes
Definition: num_vars.f90:69
num_vars::norm_disc_prec_sol
integer, public norm_disc_prec_sol
precision for normal discretization for solution
Definition: num_vars.f90:122
num_vars::ev_guess
real(dp), public ev_guess
first guess for eigenvalue
Definition: num_vars.f90:117
hdf5_vars
Variables pertaining to HDF5 and XDMF.
Definition: HDF5_vars.f90:4
grid_vars::n_r_eq
integer, public n_r_eq
nr. of normal points in equilibrium grid
Definition: grid_vars.f90:20
vmec_vars::jac_v_s
real(dp), dimension(:,:,:), allocatable, public jac_v_s
Coeff. of in cosine series (HM and FM) and norm. deriv.
Definition: VMEC_vars.f90:46
output_ops::print_ex_2d
Print 2-D output on a file.
Definition: output_ops.f90:47
num_vars::x_style
integer, public x_style
style for secondary mode numbers (1: prescribed, 2: fast)
Definition: num_vars.f90:95
grid_ops::calc_norm_range
integer function, public calc_norm_range(style, in_limits, eq_limits, X_limits, sol_limits, r_F_eq, r_F_X, r_F_sol, jq)
Calculates normal range for the input grid, the equilibrium grid and/or the solution grid.
Definition: grid_ops.f90:46
num_vars::debug_version
logical, public debug_version
debug version used
Definition: num_vars.f90:62
str_utilities
Operations on strings.
Definition: str_utilities.f90:4
num_vars::plot_magn_grid
logical, public plot_magn_grid
whether to plot the grid in real coordinates
Definition: num_vars.f90:103
vmec_vars::b_v_sub_s
real(dp), dimension(:,:,:), allocatable, public b_v_sub_s
Coeff. of B_i in cosine series (r,theta,phi) (FM)
Definition: VMEC_vars.f90:48
vmec_vars::flux_p_v
real(dp), dimension(:,:), allocatable, public flux_p_v
poloidal flux
Definition: VMEC_vars.f90:35
x_vars::min_sec_x
integer, public min_sec_x
m_X (pol. flux) or n_X (tor. flux) (only for X style 1)
Definition: X_vars.f90:127
eq_vars::max_flux_e
real(dp), public max_flux_e
max. flux in Equilibrium coordinates, set in calc_norm_range_PB3D_in
Definition: eq_vars.f90:49
str_utilities::r2strt
elemental character(len=max_str_ln) function, public r2strt(k)
Convert a real (double) to string.
Definition: str_utilities.f90:54
num_vars::prog_style
integer, public prog_style
program style (1: PB3D, 2: PB3D_POST)
Definition: num_vars.f90:53
eq_vars::r_0
real(dp), public r_0
independent normalization constant for nondimensionalization
Definition: eq_vars.f90:42
helena_ops::read_hel
integer function, public read_hel(n_r_in, use_pol_flux_H)
Reads the HELENA equilibrium data.
Definition: HELENA_ops.f90:84
num_vars::min_zeta_plot
real(dp), public min_zeta_plot
min. of zeta_plot
Definition: num_vars.f90:161
num_vars::max_it_slepc
integer, public max_it_slepc
maximum nr. of iterations for SLEPC
Definition: num_vars.f90:101
helena_vars::h_h_33
real(dp), dimension(:,:), allocatable, public h_h_33
upper metric factor (1 / gem12)
Definition: HELENA_vars.f90:32
num_vars::compare_tor_pos
logical, public compare_tor_pos
compare quantities at toroidal positions (only for POST)
Definition: num_vars.f90:151
helena_vars::r_h
real(dp), dimension(:,:), allocatable, public r_h
major radius (xout)
Definition: HELENA_vars.f90:33
grid_vars::n_r_sol
integer, public n_r_sol
nr. of normal points in solution grid
Definition: grid_vars.f90:22
num_vars::max_zvac_plot
real(dp), public max_zvac_plot
max. of Z in which to plot vacuum
Definition: num_vars.f90:168
helena_vars::flux_t_h
real(dp), dimension(:,:), allocatable, public flux_t_h
toroidal flux
Definition: HELENA_vars.f90:25
num_vars::norm_disc_style_sol
integer, public norm_disc_style_sol
style for normal discretization for solution (1: central fin. diff., 2: left fin. diff....
Definition: num_vars.f90:123
hdf5_vars::max_dim_var_1d
integer, parameter, public max_dim_var_1d
maximum dimension of var_1D
Definition: HDF5_vars.f90:21
helena_vars::h_h_12
real(dp), dimension(:,:), allocatable, public h_h_12
upper metric factor (gem12)
Definition: HELENA_vars.f90:31
vmec_vars::b_v_sub_c
real(dp), dimension(:,:,:), allocatable, public b_v_sub_c
Coeff. of B_i in sine series (r,theta,phi) (FM)
Definition: VMEC_vars.f90:47
num_vars::min_theta_plot
real(dp), public min_theta_plot
min. of theta_plot
Definition: num_vars.f90:159
grid_vars::min_par_x
real(dp), public min_par_x
min. of parallel coordinate [ ] in field-aligned grid
Definition: grid_vars.f90:24
x_vars::prim_x
integer, public prim_x
n_X (pol. flux) or m_X (tor. flux)
Definition: X_vars.f90:126
x_vars::max_r_sol
real(dp), public max_r_sol
max. normal range for pert.
Definition: X_vars.f90:136
num_vars::ev_bc
real(dp), public ev_bc
value of artificial Eigenvalue for boundary condition
Definition: num_vars.f90:116
num_vars::rich_restart_lvl
integer, public rich_restart_lvl
starting Richardson level (0: none [default])
Definition: num_vars.f90:173
helena_vars::z_h
real(dp), dimension(:,:), allocatable, public z_h
height (yout)
Definition: HELENA_vars.f90:34
hdf5_vars::var_1d_type
1D equivalent of multidimensional variables, used for internal HDF5 storage.
Definition: HDF5_vars.f90:48
output_ops::draw_ex
subroutine, public draw_ex(var_names, draw_name, nplt, draw_dim, plot_on_screen, ex_plot_style, data_name, draw_ops, extra_ops, is_animated, ranges, delay, persistent)
Use external program to draw a plot.
Definition: output_ops.f90:1079
grid_vars::n_r_in
integer, public n_r_in
nr. of normal points in input grid
Definition: grid_vars.f90:19
helena_vars::pres_h
real(dp), dimension(:,:), allocatable, public pres_h
pressure profile
Definition: HELENA_vars.f90:26
vmec_vars::flux_t_v
real(dp), dimension(:,:), allocatable, public flux_t_v
toroidal flux
Definition: VMEC_vars.f90:34
num_vars::max_nr_backtracks_hh
integer, public max_nr_backtracks_hh
maximum number of backtracks for Householder, relax. factors
Definition: num_vars.f90:132
num_vars::min_r_plot
real(dp), public min_r_plot
min. of r_plot
Definition: num_vars.f90:163
num_vars::eq_style
integer, public eq_style
either 1 (VMEC) or 2 (HELENA)
Definition: num_vars.f90:89
num_vars::plot_b
logical, public plot_b
whether to plot the magnetic field in real coordinates
Definition: num_vars.f90:104
num_vars::max_it_rich
integer, public max_it_rich
number of levels for Richardson extrapolation
Definition: num_vars.f90:127
hdf5_ops::print_hdf5_arrs
integer function, public print_hdf5_arrs(vars, PB3D_name, head_name, rich_lvl, disp_info, ind_print, remove_previous_arrs)
Prints a series of arrays, in the form of an array of pointers, to an HDF5 file.
Definition: HDF5_ops.f90:1132
vmec_vars::pres_v
real(dp), dimension(:,:), allocatable, public pres_v
pressure
Definition: VMEC_vars.f90:36
num_vars::rho_style
integer, public rho_style
style for equilibrium density profile
Definition: num_vars.f90:90
num_vars::tol_slepc
real(dp), dimension(:), allocatable, public tol_slepc
tolerance for SLEPC for different Richardson levels
Definition: num_vars.f90:118
x_vars
Variables pertaining to the perturbation quantities.
Definition: X_vars.f90:4
rich_ops::find_max_rich_lvl
integer function, public find_max_rich_lvl(lvl_rich_req, max_lvl_rich_file)
Probe to find out which Richardson levels are available.
Definition: rich_ops.f90:539
num_vars::use_pol_flux_e
logical, public use_pol_flux_e
whether poloidal flux is used in E coords.
Definition: num_vars.f90:113
num_vars::plot_j
logical, public plot_j
whether to plot the current in real coordinates
Definition: num_vars.f90:105
eq_vars::max_flux_f
real(dp), public max_flux_f
max. flux in Flux coordinates, set in calc_norm_range_PB3D_in
Definition: eq_vars.f90:50
helena_ops
Operations on HELENA variables.
Definition: HELENA_ops.f90:4
num_vars::plot_kappa
logical, public plot_kappa
whether to plot curvature
Definition: num_vars.f90:107
num_vars::tol_zero
real(dp), public tol_zero
tolerance for zeros
Definition: num_vars.f90:133
input_ops
Operations concerning giving input.
Definition: input_ops.f90:4
num_vars::norm_style
integer, public norm_style
style for normalization
Definition: num_vars.f90:92
helena_vars::h_h_11
real(dp), dimension(:,:), allocatable, public h_h_11
upper metric factor (gem11)
Definition: HELENA_vars.f90:30
eq_vars::rho_0
real(dp), public rho_0
independent normalization constant for nondimensionalization
Definition: eq_vars.f90:44
num_vars::magn_int_style
integer, public magn_int_style
style for magnetic integrals (1: trapezoidal, 2: Simpson 3/8)
Definition: num_vars.f90:124
num_vars::alpha_style
integer, public alpha_style
style for alpha (1: one field line, many turns, 2: many field lines, one turn)
Definition: num_vars.f90:100
vmec_vars::rot_t_v
real(dp), dimension(:,:), allocatable, public rot_t_v
rotational transform
Definition: VMEC_vars.f90:37
vmec_vars::b_v_s
real(dp), dimension(:,:), allocatable, public b_v_s
Coeff. of magnitude of B in cosine series (HM and FM)
Definition: VMEC_vars.f90:51
str_utilities::r2str
elemental character(len=max_str_ln) function, public r2str(k)
Convert a real (double) to string.
Definition: str_utilities.f90:42
num_vars::k_style
integer, public k_style
style for kinetic energy
Definition: num_vars.f90:93
messages::writo
subroutine, public writo(input_str, persistent, error, warning, alert)
Write output to file identified by output_i.
Definition: messages.f90:275
messages
Numerical utilities related to giving output.
Definition: messages.f90:4
vmec_vars::z_v_s
real(dp), dimension(:,:,:), allocatable, public z_v_s
Coeff. of in cosine series (FM) and norm. deriv.
Definition: VMEC_vars.f90:42
helena_vars
Variables that have to do with HELENA quantities.
Definition: HELENA_vars.f90:4
num_vars::max_rvac_plot
real(dp), public max_rvac_plot
max. of R in which to plot vacuum
Definition: num_vars.f90:166
input_ops::read_input_opts
integer function, public read_input_opts()
Reads input options from user-provided input file.
Definition: input_ops.f90:23
num_vars::plot_vmec_modes
logical, public plot_vmec_modes
plot VMEC modes
Definition: num_vars.f90:143
grid_vars::max_alpha
real(dp), public max_alpha
max. of field-line label [ ] in field-aligned grid
Definition: grid_vars.f90:27
vmec_vars::q_saf_v
real(dp), dimension(:,:), allocatable, public q_saf_v
safety factor
Definition: VMEC_vars.f90:38
num_vars::pb3d_name
character(len=max_str_ln), public pb3d_name
name of PB3D output file
Definition: num_vars.f90:139
vmec_vars::r_v_c
real(dp), dimension(:,:,:), allocatable, public r_v_c
Coeff. of in sine series (FM) and norm. deriv.
Definition: VMEC_vars.f90:39
x_vars::min_nm_x
integer, public min_nm_x
minimum for the high-n theory (debable)
Definition: X_vars.f90:130
input_ops::read_input_eq
integer function, public read_input_eq()
Reads the equilibrium input file if no Richardson restart.
Definition: input_ops.f90:1104
grid_vars
Variables pertaining to the different grids used.
Definition: grid_vars.f90:4
x_vars::max_sec_x
integer, public max_sec_x
m_X (pol. flux) or n_X (tor. flux) (only for\ c X style 1)
Definition: X_vars.f90:128
helena_vars::ias
integer, public ias
0 if top-bottom symmetric, 1 if not
Definition: HELENA_vars.f90:36
num_vars::plot_vac_pot
logical, public plot_vac_pot
whether to plot vacuum potential in POST
Definition: num_vars.f90:110
num_vars::eq_i
integer, parameter, public eq_i
file number of equilibrium file from VMEC or HELENA
Definition: num_vars.f90:183
messages::lvl_ud
subroutine, public lvl_ud(inc)
Increases/decreases lvl of output.
Definition: messages.f90:254
num_vars::max_r_plot
real(dp), public max_r_plot
max. of r_plot
Definition: num_vars.f90:164
x_vars::n_mod_x
integer, public n_mod_x
size of m_X (pol. flux) or n_X (tor. flux)
Definition: X_vars.f90:129
vmec_ops::read_vmec
integer function, public read_vmec(n_r_in, use_pol_flux_V)
Reads the VMEC equilibrium data.
Definition: VMEC_ops.f90:34
num_vars::prog_version
real(dp), parameter, public prog_version
version number
Definition: num_vars.f90:59
num_vars::plot_grid_style
integer, public plot_grid_style
style for POST plot grid (0: 3-D plots, 1: slab plots, 2: slab plots with folding,...
Definition: num_vars.f90:169
num_vars::min_rvac_plot
real(dp), public min_rvac_plot
min. of R in which to plot vacuum
Definition: num_vars.f90:165
vmec_vars
Variables that concern the output of VMEC.
Definition: VMEC_vars.f90:4
num_vars::pert_mult_factor_post
real(dp), public pert_mult_factor_post
factor with which to multiply perturbation strength for POST
Definition: num_vars.f90:176
vmec_vars::l_v_c
real(dp), dimension(:,:,:), allocatable, public l_v_c
Coeff. of in sine series (HM) and norm. deriv.
Definition: VMEC_vars.f90:43
rich_vars::min_n_par_x
integer, public min_n_par_x
min. of n_par_X (e.g. first value in Richardson loop)
Definition: rich_vars.f90:22
num_vars::use_pol_flux_f
logical, public use_pol_flux_f
whether poloidal flux is used in F coords.
Definition: num_vars.f90:114
eq_vars::psi_0
real(dp), public psi_0
derived normalization constant for nondimensionalization
Definition: eq_vars.f90:46
vmec_vars::r_v_s
real(dp), dimension(:,:,:), allocatable, public r_v_s
Coeff. of in cosine series (FM) and norm. deriv.
Definition: VMEC_vars.f90:40
rich_vars::rich_lvl
integer, public rich_lvl
current level of Richardson extrapolation
Definition: rich_vars.f90:19
vmec_vars::jac_v_c
real(dp), dimension(:,:,:), allocatable, public jac_v_c
Coeff. of in sine series (HM and FM) and norm. deriv.
Definition: VMEC_vars.f90:45
grid_vars::min_alpha
real(dp), public min_alpha
min. of field-line label [ ] in field-aligned grid
Definition: grid_vars.f90:26
num_vars::post_output_full
logical, public post_output_full
POST has output on full grids.
Definition: num_vars.f90:146
helena_vars::chi_h
real(dp), dimension(:), allocatable, public chi_h
poloidal angle
Definition: HELENA_vars.f90:23
grid_vars::n_alpha
integer, public n_alpha
nr. of field-lines
Definition: grid_vars.f90:23
input_ops::print_output_in
integer function, public print_output_in(data_name, remove_previous_arrs)
Print input quantities to an output file.
Definition: input_ops.f90:1175
output_ops
Operations concerning giving output, on the screen as well as in output files.
Definition: output_ops.f90:5
num_vars::plot_flux_q
logical, public plot_flux_q
whether to plot flux quantities in real coordinates
Definition: num_vars.f90:106
num_vars::solver_slepc_style
integer, public solver_slepc_style
style for solver (1: Krylov-Schur, 2: GD)
Definition: num_vars.f90:97
eq_vars::pres_0
real(dp), public pres_0
independent normalization constant for nondimensionalization
Definition: eq_vars.f90:43
rich_vars::req_min_n_par_x
integer, parameter, public req_min_n_par_x
required minimum n_par_X
Definition: rich_vars.f90:21
num_vars::sol_n_procs
integer, public sol_n_procs
nr. of MPI processes for solution with SLEPC
Definition: num_vars.f90:70
num_vars::max_zeta_plot
real(dp), public max_zeta_plot
max. of zeta_plot
Definition: num_vars.f90:162
num_vars::rank
integer, public rank
MPI rank.
Definition: num_vars.f90:68
num_vars::n_vac_plot
integer, dimension(2), public n_vac_plot
nr. of points in R and Z in vacuum
Definition: num_vars.f90:158
num_vars::max_tot_mem
real(dp), public max_tot_mem
maximum total memory for all processes [MB]
Definition: num_vars.f90:74
vmec_vars::l_v_s
real(dp), dimension(:,:,:), allocatable, public l_v_s
Coeff. of in cosine series (HM) and norm. deriv.
Definition: VMEC_vars.f90:44
num_vars::post_style
integer, public post_style
style for POST (1: extended grid, 2: B-aligned grid)
Definition: num_vars.f90:98
num_vars::matrix_slepc_style
integer, public matrix_slepc_style
style for matrix storage (1: sparse, 2: shell)
Definition: num_vars.f90:96
num_vars::tol_norm
real(dp), public tol_norm
tolerance for normal range (normalized to 0..1)
Definition: num_vars.f90:134
helena_vars::nchi
integer, public nchi
nr. of poloidal points
Definition: HELENA_vars.f90:35
vmec_vars::z_v_c
real(dp), dimension(:,:,:), allocatable, public z_v_c
Coeff. of in sine series (FM) and norm. deriv.
Definition: VMEC_vars.f90:41
num_vars::n_theta_plot
integer, public n_theta_plot
nr. of poloidal points in plot
Definition: num_vars.f90:156
rich_ops
Operations concerning Richardson extrapolation.
Definition: rich_ops.f90:4
num_vars::tol_rich
real(dp), public tol_rich
tolerance for Richardson extrapolation
Definition: num_vars.f90:128
num_vars::ev_style
integer, public ev_style
determines the method used for solving an EV problem
Definition: num_vars.f90:88
eq_vars::b_0
real(dp), public b_0
derived normalization constant for nondimensionalization
Definition: eq_vars.f90:45
grid_ops
Operations that have to do with the grids and different coordinate systems.
Definition: grid_ops.f90:4
num_vars::norm_disc_prec_x
integer, public norm_disc_prec_x
precision for normal discretization for perturbation
Definition: num_vars.f90:121