| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | module simple_api | ||
| 2 | !> Library interface to one SIMPLE run. | ||
| 3 | !> | ||
| 4 | !> The run configuration is the &config namelist held in module `params`; | ||
| 5 | !> set it with simple_read_config and/or by assigning `params` variables, | ||
| 6 | !> then call simple_init. The remaining calls repeat the steps of simple.x | ||
| 7 | !> (sampling, tracing, classification, output) and add single-orbit | ||
| 8 | !> diagnostics. All states are (s, theta, phi, v/v0, v_par/v) in reference | ||
| 9 | !> coordinates. The Python package pysimple wraps this module and `params`. | ||
| 10 | use, intrinsic :: iso_fortran_env, only: dp => real64 | ||
| 11 | use simple, only: tracer_t | ||
| 12 | use params, only: ntestpart, startmode, num_surf, sbeg, zstart, zend, & | ||
| 13 | times_lost, trap_par, perp_inv, iclass, class_passing, class_lost, & | ||
| 14 | ntimstep, bmod00, bmin, bmax, reallocate_arrays, reset_config, & | ||
| 15 | finalize_config, isw_field_type, integ_coords, field_input, config | ||
| 16 | use simple_profiles, only: read_profiles, reset_profiles, set_profile_table | ||
| 17 | use magfie_sub, only: TEST, CANFLUX, BOOZER | ||
| 18 | use field_boozer_chartmap, only: is_boozer_chartmap | ||
| 19 | use field_can_mod, only: integ_to_ref | ||
| 20 | use simple_main, only: validate_gpu_production_config | ||
| 21 | use simple_main, only: validate_orbit_model_config, init_run, & | ||
| 22 | init_sampling_field, sample_start, init_tracing_field, init_counters, & | ||
| 23 | trace_parallel, classify_parallel, trace_orbit, write_output, & | ||
| 24 | compute_canonical_frequencies_flat, invariants_from_state_flat, & | ||
| 25 | potato_invariants_to_simple_flat, states_from_invariants_flat, & | ||
| 26 | trace_to_cut_flat | ||
| 27 | |||
| 28 | implicit none | ||
| 29 | private | ||
| 30 | |||
| 31 | public :: simple_reset_config, simple_read_config, simple_set_config, & | ||
| 32 | simple_init, & | ||
| 33 | simple_sample, simple_sample_surface, simple_sample_volume, & | ||
| 34 | simple_set_start, simple_trace, simple_classify, simple_get_results, & | ||
| 35 | simple_get_classification, simple_write_output, simple_trace_orbit, & | ||
| 36 | simple_orbit_frequencies, simple_invariants_from_state, & | ||
| 37 | simple_states_from_invariants, simple_potato_invariants, & | ||
| 38 | simple_trace_to_cut, simple_reference_info, simple_wall_supported, & | ||
| 39 | simple_cylindrical_to_reference, simple_reference_to_cylindrical, & | ||
| 40 | simple_set_profile_table | ||
| 41 | |||
| 42 | integer, parameter :: FIELD_NONE = 0, FIELD_SAMPLING = 1, FIELD_TRACING = 2 | ||
| 43 | |||
| 44 | type(tracer_t), save :: tracer | ||
| 45 | integer, save :: active_field = FIELD_NONE | ||
| 46 | real(dp), save :: sampling_surface = -1.0_dp | ||
| 47 | !> reallocate_arrays resets the starting-line arrays (xstart, volstart) that | ||
| 48 | !> init_starting_surf fills, so they must be integrated again before sampling. | ||
| 49 | logical, save :: starting_line_stale = .true. | ||
| 50 | |||
| 51 | contains | ||
| 52 | |||
| 53 | ✗ | logical function simple_wall_supported() | |
| 54 | use stl_wall_intersection, only: stl_wall_available | ||
| 55 | |||
| 56 | ✗ | simple_wall_supported = stl_wall_available | |
| 57 | ✗ | end function simple_wall_supported | |
| 58 | |||
| 59 | ✗ | subroutine simple_set_profile_table(n, table, status) | |
| 60 | integer, intent(in) :: n | ||
| 61 | real(dp), intent(in) :: table(6, n) | ||
| 62 | integer, intent(out) :: status | ||
| 63 | |||
| 64 | ✗ | call set_profile_table(n, table, status) | |
| 65 | ✗ | end subroutine simple_set_profile_table | |
| 66 | |||
| 67 | !> Physical reference-coordinate metadata: lengths to metres, radial kind | ||
| 68 | !> 1 = VMEC s, 2 = chart rho. Normalized TEST and SPECTRE are unsupported. | ||
| 69 | ✗ | subroutine simple_reference_info(length_to_m, radial_kind, status) | |
| 70 | use reference_coordinates, only: ref_coords | ||
| 71 | use libneo_coordinates, only: chartmap_coordinate_system_t | ||
| 72 | use simple_main, only: chartmap_cart_scale_to_m | ||
| 73 | use magfie_sub, only: SPECTRE | ||
| 74 | |||
| 75 | real(dp), intent(out) :: length_to_m | ||
| 76 | integer, intent(out) :: radial_kind, status | ||
| 77 | |||
| 78 | ✗ | length_to_m = 0.01_dp | |
| 79 | ✗ | radial_kind = 1 | |
| 80 | ✗ | status = 1 | |
| 81 | ✗ | if (isw_field_type == TEST .or. isw_field_type == SPECTRE) return | |
| 82 | ✗ | if (.not. allocated(ref_coords)) return | |
| 83 | select type (cs => ref_coords) | ||
| 84 | class is (chartmap_coordinate_system_t) | ||
| 85 | ✗ | length_to_m = chartmap_cart_scale_to_m | |
| 86 | ✗ | if (length_to_m <= 0.0_dp) return | |
| 87 | ✗ | radial_kind = 2 | |
| 88 | class default | ||
| 89 | continue | ||
| 90 | end select | ||
| 91 | ✗ | status = 0 | |
| 92 | ✗ | end subroutine simple_reference_info | |
| 93 | |||
| 94 | !> Convert physical (R [m], phi [rad], Z [m]) to native reference coordinates. | ||
| 95 | !> Return one inverse status per point; invalid results are NaN. | ||
| 96 | ✗ | subroutine simple_cylindrical_to_reference(n, cylindrical, reference, status) | |
| 97 | ✗ | use, intrinsic :: ieee_arithmetic, only: ieee_is_finite, ieee_value, & | |
| 98 | ieee_quiet_nan | ||
| 99 | use reference_coordinates, only: ref_coords | ||
| 100 | |||
| 101 | integer, intent(in) :: n | ||
| 102 | real(dp), intent(in) :: cylindrical(3, n) | ||
| 103 | real(dp), intent(out) :: reference(3, n) | ||
| 104 | integer, intent(out) :: status(n) | ||
| 105 | real(dp) :: scale, point(3), result(3) | ||
| 106 | integer :: i, radial_kind, info_status | ||
| 107 | |||
| 108 | ✗ | reference = ieee_value(0.0_dp, ieee_quiet_nan) | |
| 109 | ✗ | status = -1 | |
| 110 | ✗ | call simple_reference_info(scale, radial_kind, info_status) | |
| 111 | ✗ | if (info_status /= 0) return | |
| 112 | ✗ | do i = 1, n | |
| 113 | ✗ | if (.not. all(ieee_is_finite(cylindrical(:, i)))) cycle | |
| 114 | ✗ | if (cylindrical(1, i) <= 0.0_dp) cycle | |
| 115 | ✗ | point = cylindrical(:, i) | |
| 116 | ✗ | point(1) = point(1)/scale | |
| 117 | ✗ | point(3) = point(3)/scale | |
| 118 | ✗ | call ref_coords%from_cyl(point, result, status(i)) | |
| 119 | ✗ | if (status(i) /= 0) cycle | |
| 120 | ✗ | if (.not. all(ieee_is_finite(result))) then | |
| 121 | ✗ | status(i) = -1 | |
| 122 | ✗ | cycle | |
| 123 | end if | ||
| 124 | ✗ | reference(:, i) = result | |
| 125 | end do | ||
| 126 | ✗ | end subroutine simple_cylindrical_to_reference | |
| 127 | |||
| 128 | !> Convert native reference coordinates to physical (R [m], phi [rad], Z [m]). | ||
| 129 | ✗ | subroutine simple_reference_to_cylindrical(n, reference, cylindrical, status) | |
| 130 | ✗ | use, intrinsic :: ieee_arithmetic, only: ieee_is_finite, ieee_value, & | |
| 131 | ieee_quiet_nan | ||
| 132 | use reference_coordinates, only: ref_coords | ||
| 133 | |||
| 134 | integer, intent(in) :: n | ||
| 135 | real(dp), intent(in) :: reference(3, n) | ||
| 136 | real(dp), intent(out) :: cylindrical(3, n) | ||
| 137 | integer, intent(out) :: status(n) | ||
| 138 | real(dp) :: scale, point(3) | ||
| 139 | integer :: i, radial_kind, info_status | ||
| 140 | |||
| 141 | ✗ | cylindrical = ieee_value(0.0_dp, ieee_quiet_nan) | |
| 142 | ✗ | status = -1 | |
| 143 | ✗ | call simple_reference_info(scale, radial_kind, info_status) | |
| 144 | ✗ | if (info_status /= 0) return | |
| 145 | ✗ | do i = 1, n | |
| 146 | ✗ | if (.not. all(ieee_is_finite(reference(:, i)))) cycle | |
| 147 | ✗ | call ref_coords%evaluate_cyl(reference(:, i), point) | |
| 148 | ✗ | if (.not. all(ieee_is_finite(point))) cycle | |
| 149 | ✗ | point(1) = point(1)*scale | |
| 150 | ✗ | point(3) = point(3)*scale | |
| 151 | ✗ | cylindrical(:, i) = point | |
| 152 | ✗ | status(i) = 0 | |
| 153 | end do | ||
| 154 | ✗ | end subroutine simple_reference_to_cylindrical | |
| 155 | |||
| 156 | !> Restore all &config and &profiles values to their compiled defaults. | ||
| 157 | 10 | subroutine simple_reset_config() | |
| 158 | 10 | call reset_config | |
| 159 | 10 | call reset_profiles | |
| 160 | ✗ | end subroutine simple_reset_config | |
| 161 | |||
| 162 | !> Read the &config and &profiles groups of a simple.in file. Values not | ||
| 163 | !> in the file keep their current setting. status is nonzero on error. | ||
| 164 | 2 | subroutine simple_read_config(config_file, status, message) | |
| 165 | character(*), intent(in) :: config_file | ||
| 166 | integer, intent(out) :: status | ||
| 167 | character(len=256), intent(out) :: message | ||
| 168 | |||
| 169 | integer :: unit | ||
| 170 | |||
| 171 | 2 | message = '' | |
| 172 | open (newunit=unit, file=config_file, status='old', action='read', & | ||
| 173 | 2 | iostat=status, iomsg=message) | |
| 174 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
2 | if (status /= 0) return |
| 175 | 2 | read (unit, nml=config, iostat=status, iomsg=message) | |
| 176 | 2 | close (unit) | |
| 177 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
2 | if (status /= 0) return |
| 178 | 2 | call read_profiles(config_file) | |
| 179 | end subroutine simple_read_config | ||
| 180 | |||
| 181 | !> Set &config values from namelist text such as "&config ntestpart=8 /". | ||
| 182 | !> status is nonzero for unknown names or malformed values. | ||
| 183 | ✗ | subroutine simple_set_config(namelist_text, status, message) | |
| 184 | character(*), intent(in) :: namelist_text | ||
| 185 | integer, intent(out) :: status | ||
| 186 | character(len=256), intent(out) :: message | ||
| 187 | |||
| 188 | ✗ | message = '' | |
| 189 | ✗ | read (namelist_text, nml=config, iostat=status, iomsg=message) | |
| 190 | ✗ | end subroutine simple_set_config | |
| 191 | |||
| 192 | !> Validate the configuration and load the field, as simple.x does before | ||
| 193 | !> sampling. | ||
| 194 | !> A Boozer chartmap input selects isw_field_type = 2 (Boozer) when neither | ||
| 195 | !> isw_field_type nor integ_coords was set. | ||
| 196 | 12 | subroutine simple_init() | |
| 197 | logical :: field_type_unset | ||
| 198 | |||
| 199 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
|
12 | field_type_unset = integ_coords == -1000 .and. isw_field_type == CANFLUX |
| 200 | 12 | call finalize_config | |
| 201 | ✗ | if (field_type_unset .and. len_trim(field_input) > 0) then | |
| 202 | ✗ | if (is_boozer_chartmap(field_input)) then | |
| 203 | ✗ | isw_field_type = BOOZER | |
| 204 | ✗ | integ_coords = BOOZER | |
| 205 | end if | ||
| 206 | end if | ||
| 207 | 12 | call validate_orbit_model_config | |
| 208 | 12 | call validate_gpu_production_config | |
| 209 | 12 | call init_run(tracer) | |
| 210 | 12 | call init_sampling_field | |
| 211 | 12 | starting_line_stale = .false. | |
| 212 | 12 | active_field = FIELD_SAMPLING | |
| 213 | 12 | sampling_surface = sbeg(1) | |
| 214 | 12 | end subroutine simple_init | |
| 215 | |||
| 216 | !> Sample n start states as configured by startmode, sbeg and num_surf. | ||
| 217 | 13 | subroutine simple_sample(n, z) | |
| 218 | integer, intent(in) :: n | ||
| 219 | real(dp), intent(out) :: z(5, n) | ||
| 220 | |||
| 221 | 13 | call allocate_particles(n) | |
| 222 | 13 | call use_sampling_field(needs_starting_line=.true.) | |
| 223 | 13 | call sample_start | |
| 224 |
4/4✓ Branch 0 taken 700040 times.
✓ Branch 1 taken 12 times.
✓ Branch 2 taken 3500200 times.
✓ Branch 3 taken 700040 times.
|
4200252 | z = zstart(1:5, 1:n) |
| 225 | 12 | end subroutine simple_sample | |
| 226 | |||
| 227 | !> Sample n start states on the flux surface s along the starting field line. | ||
| 228 | 2 | subroutine simple_sample_surface(n, s, z) | |
| 229 | integer, intent(in) :: n | ||
| 230 | real(dp), intent(in) :: s | ||
| 231 | real(dp), intent(out) :: z(5, n) | ||
| 232 | |||
| 233 | 2 | startmode = 1 | |
| 234 | 2 | num_surf = 1 | |
| 235 | 2 | sbeg(1) = s | |
| 236 | 2 | call simple_sample(n, z) | |
| 237 | 2 | end subroutine simple_sample_surface | |
| 238 | |||
| 239 | !> Sample n start states uniformly in the volume s_inner < s < s_outer. | ||
| 240 | ✗ | subroutine simple_sample_volume(n, s_inner, s_outer, z) | |
| 241 | integer, intent(in) :: n | ||
| 242 | real(dp), intent(in) :: s_inner, s_outer | ||
| 243 | real(dp), intent(out) :: z(5, n) | ||
| 244 | |||
| 245 | ✗ | startmode = 5 | |
| 246 | ✗ | num_surf = 2 | |
| 247 | ✗ | sbeg(1:2) = [s_inner, s_outer] | |
| 248 | ✗ | call simple_sample(n, z) | |
| 249 | ✗ | end subroutine simple_sample_volume | |
| 250 | |||
| 251 | !> Use the n states z as start states. When num_surf = 1 and all states | ||
| 252 | !> share one radius, that radius becomes the starting surface, so bmin and | ||
| 253 | !> bmax refer to it. | ||
| 254 | ✗ | subroutine simple_set_start(n, z) | |
| 255 | integer, intent(in) :: n | ||
| 256 | real(dp), intent(in) :: z(5, n) | ||
| 257 | |||
| 258 | ✗ | call allocate_particles(n) | |
| 259 | ✗ | zstart(1:5, 1:n) = z | |
| 260 | ✗ | if (n > 0) then | |
| 261 | ✗ | if (maxval(z(1, :)) <= minval(z(1, :))) call set_start_surface(z(1, 1)) | |
| 262 | end if | ||
| 263 | ✗ | end subroutine simple_set_start | |
| 264 | |||
| 265 | !> Trace the start states to trace_time or loss, as simple.x does. | ||
| 266 | 1 | subroutine simple_trace() | |
| 267 | 1 | call use_tracing_field | |
| 268 | 1 | call init_counters | |
| 269 | 1 | call trace_parallel(tracer) | |
| 270 | 1 | end subroutine simple_trace | |
| 271 | |||
| 272 | !> Trace the start states with orbit classification enabled. | ||
| 273 | ✗ | subroutine simple_classify() | |
| 274 | ✗ | call use_tracing_field | |
| 275 | ✗ | call init_counters | |
| 276 | ✗ | call classify_parallel(tracer) | |
| 277 | ✗ | end subroutine simple_classify | |
| 278 | |||
| 279 | !> Copy the per-particle results of the last trace or classification. | ||
| 280 | !> times_lost is trace_time for confined, -1 for skipped and NaN for | ||
| 281 | !> numerically failed particles. | ||
| 282 | ✗ | subroutine simple_get_results(n, z_end, loss_time, trap_parameter, & | |
| 283 | ✗ | perp_invariant) | |
| 284 | integer, intent(in) :: n | ||
| 285 | real(dp), intent(out) :: z_end(5, n), loss_time(n), trap_parameter(n), & | ||
| 286 | perp_invariant(n) | ||
| 287 | |||
| 288 | ✗ | z_end = zend(1:5, 1:n) | |
| 289 | ✗ | loss_time = times_lost(1:n) | |
| 290 | ✗ | trap_parameter = trap_par(1:n) | |
| 291 | ✗ | perp_invariant = perp_inv(1:n) | |
| 292 | ✗ | end subroutine simple_get_results | |
| 293 | |||
| 294 | !> Copy the classification of the last simple_classify call. passing and | ||
| 295 | !> lost are 0 or 1; classes holds the J_parallel, topological and fractal | ||
| 296 | !> class (0 unclassified, 1 regular, 2 chaotic) per particle. | ||
| 297 | ✗ | subroutine simple_get_classification(n, passing, lost, classes) | |
| 298 | integer, intent(in) :: n | ||
| 299 | integer, intent(out) :: passing(n), lost(n), classes(3, n) | ||
| 300 | |||
| 301 | ✗ | passing = merge(1, 0, class_passing(1:n)) | |
| 302 | ✗ | lost = merge(1, 0, class_lost(1:n)) | |
| 303 | ✗ | classes = iclass(1:3, 1:n) | |
| 304 | ✗ | end subroutine simple_get_classification | |
| 305 | |||
| 306 | !> Write the simple.x output files for the last trace. | ||
| 307 | 1 | subroutine simple_write_output() | |
| 308 | 1 | call write_output | |
| 309 | 1 | end subroutine simple_write_output | |
| 310 | |||
| 311 | !> Trace one orbit and return its ntimstep samples in reference coordinates. | ||
| 312 | ✗ | subroutine simple_trace_orbit(z0, nt, trajectory, times, z_end, loss_time) | |
| 313 | real(dp), intent(in) :: z0(5) | ||
| 314 | integer, intent(in) :: nt | ||
| 315 | real(dp), intent(out) :: trajectory(5, nt), times(nt), z_end(5) | ||
| 316 | real(dp), intent(out) :: loss_time | ||
| 317 | |||
| 318 | real(dp), allocatable :: integ_trajectory(:, :) | ||
| 319 | integer :: it | ||
| 320 | |||
| 321 | ✗ | if (nt /= ntimstep) error stop 'simple_trace_orbit: nt must equal ntimstep' | |
| 322 | ✗ | call simple_set_start(1, reshape(z0, [5, 1])) | |
| 323 | ✗ | call use_tracing_field | |
| 324 | ✗ | call init_counters | |
| 325 | ✗ | allocate (integ_trajectory(5, nt)) | |
| 326 | ✗ | call trace_orbit(tracer, 1, integ_trajectory, times) | |
| 327 | ✗ | do it = 1, nt | |
| 328 | ✗ | call integ_to_ref(integ_trajectory(1:3, it), trajectory(1:3, it)) | |
| 329 | ✗ | trajectory(4:5, it) = integ_trajectory(4:5, it) | |
| 330 | end do | ||
| 331 | ✗ | z_end = zend(1:5, 1) | |
| 332 | ✗ | loss_time = times_lost(1) | |
| 333 | ✗ | end subroutine simple_trace_orbit | |
| 334 | |||
| 335 | !> Bounce/transit period, toroidal shift and frequencies of one orbit. | ||
| 336 | !> values: period, period_std, delta_phi, delta_phi_std, omega_b, omega_phi. | ||
| 337 | !> metadata: status, orbit_class, parallel_direction, n_periods, n_steps. | ||
| 338 | ✗ | subroutine simple_orbit_frequencies(z0, n_periods, max_steps, values, metadata) | |
| 339 | real(dp), intent(in) :: z0(5) | ||
| 340 | integer, intent(in) :: n_periods, max_steps | ||
| 341 | real(dp), intent(out) :: values(6) | ||
| 342 | integer, intent(out) :: metadata(5) | ||
| 343 | |||
| 344 | ✗ | call set_start_surface(z0(1)) | |
| 345 | ✗ | call use_tracing_field | |
| 346 | call compute_canonical_frequencies_flat(tracer, z0, n_periods, max_steps, & | ||
| 347 | ✗ | values, metadata) | |
| 348 | ✗ | end subroutine simple_orbit_frequencies | |
| 349 | |||
| 350 | !> Invariants of one state. values: h0, j_perp, p_phi, axis flux, outward | ||
| 351 | !> sign, edge flux. | ||
| 352 | ✗ | subroutine simple_invariants_from_state(z0, values, status) | |
| 353 | real(dp), intent(in) :: z0(5) | ||
| 354 | real(dp), intent(out) :: values(6) | ||
| 355 | integer, intent(out) :: status | ||
| 356 | |||
| 357 | ✗ | call use_tracing_field | |
| 358 | ✗ | call invariants_from_state_flat(tracer, z0, values, status) | |
| 359 | ✗ | end subroutine simple_invariants_from_state | |
| 360 | |||
| 361 | !> Start states on the axisymmetric cut matching (h0, j_perp, p_phi), | ||
| 362 | !> ordered by major radius. metadata: sigma, section branch, section kind. | ||
| 363 | ✗ | subroutine simple_states_from_invariants(invariants, max_solutions, states, & | |
| 364 | ✗ | metadata, residuals, cylindrical, n_solutions, status) | |
| 365 | real(dp), intent(in) :: invariants(3) | ||
| 366 | integer, intent(in) :: max_solutions | ||
| 367 | real(dp), intent(out) :: states(5, max_solutions) | ||
| 368 | integer, intent(out) :: metadata(3, max_solutions) | ||
| 369 | real(dp), intent(out) :: residuals(max_solutions) | ||
| 370 | real(dp), intent(out) :: cylindrical(3, max_solutions) | ||
| 371 | integer, intent(out) :: n_solutions, status | ||
| 372 | |||
| 373 | ✗ | call use_tracing_field | |
| 374 | call states_from_invariants_flat(tracer, invariants, max_solutions, & | ||
| 375 | ✗ | states, metadata, residuals, cylindrical, n_solutions, status) | |
| 376 | ✗ | end subroutine simple_states_from_invariants | |
| 377 | |||
| 378 | !> Convert POTATO invariants (h0, j_perp, psi_star) to SIMPLE's convention. | ||
| 379 | ✗ | subroutine simple_potato_invariants(potato, psi_axis, psi_edge, v0_ratio, & | |
| 380 | invariants) | ||
| 381 | real(dp), intent(in) :: potato(3), psi_axis, psi_edge, v0_ratio | ||
| 382 | real(dp), intent(out) :: invariants(3) | ||
| 383 | |||
| 384 | call potato_invariants_to_simple_flat(potato, psi_axis, psi_edge, & | ||
| 385 | ✗ | v0_ratio, invariants) | |
| 386 | ✗ | end subroutine simple_potato_invariants | |
| 387 | |||
| 388 | !> Trace one orbit to its next tip (cut 0), toroidal-period cut (1) or | ||
| 389 | !> either (-1). cut_state holds the state and the parallel invariant. | ||
| 390 | ✗ | subroutine simple_trace_to_cut(z0, cut, max_events, cut_state, cut_type, & | |
| 391 | status) | ||
| 392 | real(dp), intent(in) :: z0(5) | ||
| 393 | integer, intent(in) :: cut, max_events | ||
| 394 | real(dp), intent(out) :: cut_state(6) | ||
| 395 | integer, intent(out) :: cut_type, status | ||
| 396 | |||
| 397 | ✗ | call set_start_surface(z0(1)) | |
| 398 | ✗ | call use_tracing_field | |
| 399 | call trace_to_cut_flat(tracer, z0, cut, max_events, cut_state, cut_type, & | ||
| 400 | ✗ | status) | |
| 401 | ✗ | end subroutine simple_trace_to_cut | |
| 402 | |||
| 403 | 13 | subroutine allocate_particles(n) | |
| 404 | integer, intent(in) :: n | ||
| 405 | |||
| 406 | 13 | ntestpart = n | |
| 407 | 13 | call reallocate_arrays | |
| 408 | 13 | starting_line_stale = .true. | |
| 409 | end subroutine allocate_particles | ||
| 410 | |||
| 411 | ✗ | subroutine set_start_surface(s) | |
| 412 | real(dp), intent(in) :: s | ||
| 413 | |||
| 414 | ✗ | if (isw_field_type == TEST) then | |
| 415 | ! Analytic tokamak: B = 1 - r cos(theta), so B ranges over 1 -+ r. | ||
| 416 | ✗ | bmod00 = 1.0_dp | |
| 417 | ✗ | bmin = 1.0_dp - s | |
| 418 | ✗ | bmax = 1.0_dp + s | |
| 419 | ✗ | else if (num_surf == 1) then | |
| 420 | ✗ | sbeg(1) = s | |
| 421 | ✗ | call use_sampling_field | |
| 422 | end if | ||
| 423 | ✗ | end subroutine set_start_surface | |
| 424 | |||
| 425 | 13 | subroutine use_sampling_field(needs_starting_line) | |
| 426 | ! The starting field line depends on sbeg(1); integrate it only when the | ||
| 427 | ! sampling field is not active on that surface, or when the caller | ||
| 428 | ! samples and the line arrays were reset by reallocate_arrays. | ||
| 429 | logical, intent(in), optional :: needs_starting_line | ||
| 430 | logical :: refill | ||
| 431 | |||
| 432 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 13 times.
|
13 | if (active_field == FIELD_NONE) error stop 'simple_api: call simple_init first' |
| 433 | 13 | refill = .false. | |
| 434 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 13 times.
|
13 | if (present(needs_starting_line)) refill = needs_starting_line .and. & |
| 435 |
2/4✗ Branch 0 not taken.
✓ Branch 1 taken 13 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 13 times.
|
13 | starting_line_stale |
| 436 |
2/4✓ Branch 0 taken 13 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 13 times.
|
13 | if (active_field == FIELD_SAMPLING .and. .not. refill) then |
| 437 | ✗ | if (abs(sbeg(1) - sampling_surface) <= 0.0_dp) return | |
| 438 | end if | ||
| 439 | 13 | call init_sampling_field | |
| 440 | 13 | starting_line_stale = .false. | |
| 441 | 13 | active_field = FIELD_SAMPLING | |
| 442 | 13 | sampling_surface = sbeg(1) | |
| 443 | end subroutine use_sampling_field | ||
| 444 | |||
| 445 | 1 | subroutine use_tracing_field | |
| 446 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 1 time.
|
1 | if (active_field == FIELD_NONE) error stop 'simple_api: call simple_init first' |
| 447 |
1/2✓ Branch 0 taken 1 time.
✗ Branch 1 not taken.
|
1 | if (active_field /= FIELD_TRACING) call init_tracing_field |
| 448 | 1 | active_field = FIELD_TRACING | |
| 449 | 1 | end subroutine use_tracing_field | |
| 450 | |||
| 451 | ✗ | end module simple_api | |
| 452 |