| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | ✗ | program diag_orbit_main | |
| 2 | !> Diagnostic application for orbit trajectory analysis | ||
| 3 | !> | ||
| 4 | !> Reads configuration from simple.in (default) or specified file, | ||
| 5 | !> validates the integration setup, initializes the field exactly like simple.x, | ||
| 6 | !> and provides detailed orbit trajectory plotting for the Nth particle | ||
| 7 | |||
| 8 | use params, only: read_config, netcdffile, ns_s, ns_tp, multharm, & | ||
| 9 | integmode, params_init, isw_field_type, dtaumin, relerr, ntestpart, ntimstep, ntau, & | ||
| 10 | ✗ | zstart, startmode, grid_density, special_ants_file, reuse_batch, num_surf, sbeg | |
| 11 | use simple, only: tracer_t | ||
| 12 | use simple_main, only: init_field | ||
| 13 | use magfie_sub, only: init_magfie, VMEC | ||
| 14 | use samplers, only: init_starting_surf, sample, START_FILE | ||
| 15 | use timing, only: init_timer, print_phase_time | ||
| 16 | use diag_orbit, only: integrate_orbit_with_trajectory_debug | ||
| 17 | use orbit_symplectic_base, only: symplectic_integrator_t | ||
| 18 | use field_can_mod, only: field_can_t | ||
| 19 | |||
| 20 | implicit none | ||
| 21 | |||
| 22 | character(256) :: config_file, particle_arg | ||
| 23 | type(tracer_t) :: norb | ||
| 24 | type(symplectic_integrator_t) :: si | ||
| 25 | type(field_can_t) :: field_can | ||
| 26 | integer :: particle_number | ||
| 27 | |||
| 28 | ! Initialize timing | ||
| 29 | ✗ | call init_timer() | |
| 30 | |||
| 31 | ! Read configuration file name from command line arguments | ||
| 32 | ✗ | if (command_argument_count() == 0) then | |
| 33 | ✗ | config_file = 'simple.in' | |
| 34 | ✗ | particle_number = 1 ! Default to first particle | |
| 35 | ✗ | elseif (command_argument_count() == 1) then | |
| 36 | ✗ | call get_command_argument(1, particle_arg) | |
| 37 | ✗ | read (particle_arg, *) particle_number | |
| 38 | ✗ | config_file = 'simple.in' | |
| 39 | ✗ | elseif (command_argument_count() == 2) then | |
| 40 | ✗ | call get_command_argument(1, config_file) | |
| 41 | ✗ | call get_command_argument(2, particle_arg) | |
| 42 | ✗ | read (particle_arg, *) particle_number | |
| 43 | else | ||
| 44 | ✗ | print *, 'Usage: ./diag_orbit.x [config_file] [particle_number]' | |
| 45 | ✗ | print *, ' or: ./diag_orbit.x [particle_number]' | |
| 46 | ✗ | print *, 'Example: ./diag_orbit.x simple.in 2' | |
| 47 | ✗ | print *, ' ./diag_orbit.x 3' | |
| 48 | ✗ | stop | |
| 49 | end if | ||
| 50 | ✗ | call print_phase_time('Command line parsing completed') | |
| 51 | |||
| 52 | ! Initialize the system following simple.x pattern | ||
| 53 | ✗ | call read_config(config_file) | |
| 54 | ✗ | call print_phase_time('Configuration reading completed') | |
| 55 | |||
| 56 | ! Validate particle number against ntestpart | ||
| 57 | ✗ | if (particle_number < 1 .or. particle_number > ntestpart) then | |
| 58 | ✗ | print *, 'ERROR: Invalid particle number!' | |
| 59 | ✗ | print '(A,I0)', 'Requested particle: ', particle_number | |
| 60 | ✗ | print '(A,I0)', 'Available particles (ntestpart): ', ntestpart | |
| 61 | ✗ | print *, 'Please adjust particle number or ntestpart in config file.' | |
| 62 | ✗ | error stop 'Invalid particle number for orbit trajectory diagnostic' | |
| 63 | end if | ||
| 64 | ✗ | call print_phase_time('Particle number validation completed') | |
| 65 | |||
| 66 | ! Use the complete field initialization from simple_main | ||
| 67 | ✗ | call init_field(norb, netcdffile, ns_s, ns_tp, multharm, integmode) | |
| 68 | ✗ | call print_phase_time('Complete field initialization completed') | |
| 69 | |||
| 70 | ✗ | call params_init | |
| 71 | ✗ | call print_phase_time('Parameter initialization completed') | |
| 72 | |||
| 73 | ! Initialize VMEC magnetic field (required for sampling) | ||
| 74 | ✗ | call init_magfie(VMEC) | |
| 75 | ✗ | call print_phase_time('VMEC magnetic field initialization completed') | |
| 76 | |||
| 77 | ! Initialize starting surfaces (required before sampling) | ||
| 78 | ✗ | call init_starting_surf | |
| 79 | ✗ | call print_phase_time('Starting surface initialization completed') | |
| 80 | |||
| 81 | ! Perform particle sampling exactly like simple_main.f90 | ||
| 82 | ✗ | if (1 == startmode) then | |
| 83 | ✗ | if ((0d0 < grid_density) .and. (1d0 > grid_density)) then | |
| 84 | ✗ | call sample(zstart, grid_density) | |
| 85 | else | ||
| 86 | ✗ | call sample(zstart) | |
| 87 | end if | ||
| 88 | ✗ | elseif (2 == startmode) then | |
| 89 | ✗ | call sample(zstart, START_FILE) | |
| 90 | ✗ | elseif (3 == startmode) then | |
| 91 | ✗ | call sample(special_ants_file) | |
| 92 | ✗ | elseif (4 == startmode) then | |
| 93 | ✗ | call sample(zstart, reuse_batch) | |
| 94 | ✗ | elseif (5 == startmode) then | |
| 95 | ✗ | if (0 == num_surf) then | |
| 96 | ✗ | call sample(zstart, 0.0d0, 1.0d0) | |
| 97 | ✗ | elseif (1 == num_surf) then | |
| 98 | ✗ | call sample(zstart, 0.0d0, sbeg(1)) | |
| 99 | ✗ | elseif (2 == num_surf) then | |
| 100 | ✗ | call sample(zstart, sbeg(1), sbeg(num_surf)) | |
| 101 | else | ||
| 102 | ✗ | print *, 'Invalid surface range for volume sample defined.' | |
| 103 | ✗ | error stop 'Invalid surface range for volume sample' | |
| 104 | end if | ||
| 105 | else | ||
| 106 | ✗ | print *, 'Invalid startmode: ', startmode | |
| 107 | ✗ | error stop 'Invalid startmode' | |
| 108 | end if | ||
| 109 | ✗ | call print_phase_time('Particle sampling completed') | |
| 110 | |||
| 111 | ✗ | print *, "Orbit Trajectory Diagnostic Program" | |
| 112 | ✗ | print *, "===================================" | |
| 113 | ✗ | print * | |
| 114 | ✗ | print *, "This program provides detailed orbit trajectory plotting" | |
| 115 | ✗ | print *, "for individual particles using real physics integration." | |
| 116 | ✗ | print * | |
| 117 | ✗ | print '(A,A)', "Configuration file: ", trim(config_file) | |
| 118 | ✗ | print '(A,I0,A,I0)', "Selected particle: ", particle_number, " out of ", ntestpart | |
| 119 | ✗ | print '(A,I0)', "Field type (isw_field_type): ", isw_field_type | |
| 120 | ✗ | print '(A,I0)', "Integration mode: ", integmode | |
| 121 | ✗ | print '(A,I0,A,I0,A,I0,A)', "Integration: ", ntimstep, " macrosteps × ", ntau, " substeps = ", ntimstep*ntau, " total timesteps" | |
| 122 | ✗ | print * | |
| 123 | |||
| 124 | ✗ | print '(A,ES12.5)', 'dtaumin (integration time step): ', dtaumin | |
| 125 | ✗ | print '(A,I0)', 'ntau (substeps per dtau): ', ntau | |
| 126 | ✗ | print '(A,ES12.5)', 'Relative tolerance: ', relerr | |
| 127 | ✗ | print * | |
| 128 | |||
| 129 | ! NOTE: Symplectic integrator initialization happens per-particle in integrate_orbit_with_trajectory_debug | ||
| 130 | ! following the exact same sequence as simple_main.f90 trace_orbit() | ||
| 131 | |||
| 132 | ! Perform orbit trajectory diagnostic integration | ||
| 133 | ✗ | call integrate_orbit_with_trajectory_debug(si, field_can, particle_number) | |
| 134 | |||
| 135 | ✗ | call print_phase_time('Orbit trajectory diagnostic analysis completed') | |
| 136 | |||
| 137 | ✗ | print * | |
| 138 | ✗ | print *, "Orbit Trajectory Diagnostic completed successfully!" | |
| 139 | ✗ | print *, "Generated detailed trajectory plots showing real particle" | |
| 140 | ✗ | print *, "motion in the magnetic field using symplectic integration." | |
| 141 | |||
| 142 | ✗ | end program diag_orbit_main | |
| 143 |