| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | module field_newton | ||
| 2 | !> Module for field-agnostic Newton iterations | ||
| 3 | !> Provides Newton solvers for coordinate transformations | ||
| 4 | |||
| 5 | use, intrinsic :: iso_fortran_env, only : dp => real64 | ||
| 6 | use field_base, only : magnetic_field_t | ||
| 7 | |||
| 8 | implicit none | ||
| 9 | private | ||
| 10 | |||
| 11 | public :: newton_theta_from_canonical | ||
| 12 | |||
| 13 | contains | ||
| 14 | |||
| 15 | !> Newton iteration to find field-specific theta from canonical theta | ||
| 16 | !> Solves: theta + Lambda(s, theta, phi) = vartheta_canonical | ||
| 17 | !> where Lambda is the stream function that depends on the field type | ||
| 18 | 44257846 | subroutine newton_theta_from_canonical(field, s, vartheta_canonical, varphi, & | |
| 19 | theta, converged, iter_count) | ||
| 20 | class(magnetic_field_t), intent(in) :: field | ||
| 21 | real(dp), intent(in) :: s ! Normalized flux coordinate | ||
| 22 | real(dp), intent(in) :: vartheta_canonical ! Canonical poloidal angle | ||
| 23 | real(dp), intent(in) :: varphi ! Toroidal angle | ||
| 24 | real(dp), intent(inout) :: theta ! Field-specific theta (initial guess on input) | ||
| 25 | logical, intent(out) :: converged ! Convergence flag | ||
| 26 | integer, intent(out), optional :: iter_count ! Number of iterations | ||
| 27 | |||
| 28 | ! Local variables | ||
| 29 | real(dp), parameter :: epserr = 1.0e-12_dp ! Convergence tolerance | ||
| 30 | integer, parameter :: max_iter = 100 | ||
| 31 | real(dp) :: alam, dl_dt, deltheta | ||
| 32 | integer :: iter | ||
| 33 | |||
| 34 | ! Initialize | ||
| 35 | 9403985 | converged = .false. | |
| 36 | |||
| 37 | ! Newton iteration loop | ||
| 38 |
1/2✓ Branch 0 taken 34853861 times.
✗ Branch 1 not taken.
|
34853861 | do iter = 1, max_iter |
| 39 | ! Get stream function and its theta derivative for current field | ||
| 40 | 34853861 | call get_stream_function(field, s, theta, varphi, alam, dl_dt) | |
| 41 | |||
| 42 | ! Newton update: solve vartheta = theta + Lambda | ||
| 43 | ! f(theta) = theta + Lambda(theta) - vartheta = 0 | ||
| 44 | ! f'(theta) = 1 + dLambda/dtheta | ||
| 45 | 34853861 | deltheta = (vartheta_canonical - theta - alam) / (1.0_dp + dl_dt) | |
| 46 | 34853861 | theta = theta + deltheta | |
| 47 | |||
| 48 | ! Check convergence | ||
| 49 |
2/2✓ Branch 1 taken 9403985 times.
✓ Branch 2 taken 25449876 times.
|
34853861 | if (abs(deltheta) < epserr) then |
| 50 | 9403985 | converged = .true. | |
| 51 | 9403985 | exit | |
| 52 | end if | ||
| 53 | end do | ||
| 54 | |||
| 55 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 9403985 times.
|
9403985 | if (present(iter_count)) iter_count = iter |
| 56 | |||
| 57 | ! Warn if not converged | ||
| 58 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 9403985 times.
|
9403985 | if (.not. converged) then |
| 59 | ✗ | print *, 'WARNING: Newton iteration for theta did not converge' | |
| 60 | ✗ | print *, ' s =', s, ', vartheta =', vartheta_canonical, ', varphi =', varphi | |
| 61 | ✗ | print *, ' Final theta =', theta, ', deltheta =', deltheta | |
| 62 | end if | ||
| 63 | |||
| 64 | 9403985 | end subroutine newton_theta_from_canonical | |
| 65 | |||
| 66 | |||
| 67 | !> Get stream function Lambda and its derivative for a given field | ||
| 68 | !> This abstracts the field-specific stream function evaluation | ||
| 69 | 34853861 | subroutine get_stream_function(mag_field, s, theta, varphi, alam, dl_dt) | |
| 70 | use field, only: vmec_field_t | ||
| 71 | use vmec_field_eval, only: vmec_lambda_interpolate_with_field | ||
| 72 | |||
| 73 | class(magnetic_field_t), intent(in) :: mag_field | ||
| 74 | real(dp), intent(in) :: s, theta, varphi | ||
| 75 | real(dp), intent(out) :: alam ! Stream function Lambda | ||
| 76 | real(dp), intent(out) :: dl_dt ! dLambda/dtheta | ||
| 77 | |||
| 78 | ! Dispatch based on field type | ||
| 79 | select type (mag_field) | ||
| 80 | type is (vmec_field_t) | ||
| 81 | ! For VMEC, use the existing lambda interpolation | ||
| 82 | 34853861 | call vmec_lambda_interpolate_with_field(mag_field, s, theta, varphi, alam, dl_dt) | |
| 83 | |||
| 84 | class default | ||
| 85 | ! For other fields, stream function may not be defined | ||
| 86 | ✗ | alam = 0.0_dp | |
| 87 | ✗ | dl_dt = 0.0_dp | |
| 88 | end select | ||
| 89 | |||
| 90 | 34853861 | end subroutine get_stream_function | |
| 91 | |||
| 92 |
1/2✓ Branch 0 taken 34853861 times.
✗ Branch 1 not taken.
|
34853861 | end module field_newton |
| 93 |