forked from landreman/regcoil
-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathcompute_offset_surface_mod.f90
78 lines (55 loc) · 2.68 KB
/
compute_offset_surface_mod.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
module compute_offset_surface_mod
use stel_kinds
implicit none
private
public :: compute_offset_surface_xyz_of_thetazeta
real(dp) :: theta_rootSolve, zeta_rootSolve_target, separation
!$OMP threadprivate(theta_rootSolve, zeta_rootSolve_target, separation)
contains
subroutine compute_offset_surface_xyz_of_thetazeta(theta_rootSolve_in,zeta_rootSolve_target_in,x_offsetSurface,y_offsetSurface,z_offsetSurface,separation_in)
use stel_kinds
implicit none
real(dp), intent(in) :: theta_rootSolve_in, zeta_rootSolve_target_in, separation_in
real(dp), intent(out) :: x_offsetSurface, y_offsetSurface, z_offsetSurface
integer :: fzeroFlag
real(dp) :: rootSolve_abserr, rootSolve_relerr, zeta_rootSolve_min, zeta_rootSolve_max
real(dp) :: zeta_plasma_rootSolveSolution
theta_rootSolve = theta_rootSolve_in
zeta_rootSolve_target = zeta_rootSolve_target_in
separation = separation_in
!rootSolve_abserr = 0
!rootSolve_relerr = 0
rootSolve_abserr = 1.0e-10_dp
rootSolve_relerr = 1.0e-10_dp
zeta_rootSolve_min = zeta_rootSolve_target - 0.3
zeta_rootSolve_max = zeta_rootSolve_target + 0.3
call fzero(fzero_residual, zeta_rootSolve_min, zeta_rootSolve_max, zeta_rootSolve_target, &
rootSolve_relerr, rootSolve_abserr, fzeroFlag)
! Note: fzero returns its answer in zeta_rootSolve_min
zeta_plasma_rootSolveSolution = zeta_rootSolve_min
if (fzeroFlag == 4) then
stop "ERROR: fzero returned error 4: no sign change in residual"
else if (fzeroFlag > 2) then
print *,"WARNING: fzero returned an error code:",fzeroFlag
end if
call expand_plasma_surface(theta_rootSolve, zeta_plasma_rootSolveSolution, separation, x_offsetSurface, y_offsetSurface, z_offsetSurface)
end subroutine compute_offset_surface_xyz_of_thetazeta
!------------------------------------------------------------------------------------
function fzero_residual(zeta_plasma_test)
use global_variables, only: nfp
use stel_constants
implicit none
real(dp) :: zeta_plasma_test, fzero_residual
real(dp) :: x_outer, y_outer, z_outer, zeta_outer_new, zeta_error
call expand_plasma_surface(theta_rootSolve, zeta_plasma_test, separation, x_outer, y_outer, z_outer)
zeta_outer_new = atan2(y_outer,x_outer)
zeta_error = zeta_outer_new - zeta_rootSolve_target
if (zeta_error < -pi) then
zeta_error = zeta_error + twopi
end if
if (zeta_error > pi) then
zeta_error = zeta_error - twopi
end if
fzero_residual = zeta_error
end function fzero_residual
end module compute_offset_surface_mod