on one flux surface, find for a list of in thet_j,zeta_j, the corresponding (thet_j,zeta_j) positions, given Here, new PEST angles are theta=theta+lambda(theta,zeta) zeta=zeta, so a 1D root search in theta is is enough
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| integer, | intent(in) | :: | nrho |
number of surfaces, (second dimension of LA_in and nu_in modes) |
||
| type(t_fBase), | intent(in) | :: | fbase_in | |||
| real(kind=wp), | intent(in) | :: | LA_in(1:fbase_in%modes,nrho) | |||
| integer, | intent(in) | :: | tz_dim | |||
| real(kind=wp), | intent(in) | :: | tz_pest(2,tz_dim) | |||
| real(kind=wp), | intent(out) | :: | thetzeta_out(2,tz_dim,nrho) |
theta,zeta position in original angles, for given pest angles |
SUBROUTINE find_pest_angles(nrho,fbase_in,LA_in,tz_dim,tz_pest,thetzeta_out) ! MODULES USE MODgvec_Globals,ONLY: UNIT_stdOut,ProgressBar,testlevel USE MODgvec_fbase ,ONLY: t_fbase IMPLICIT NONE !----------------------------------------------------------------------------------------------------------------------------------- ! INPUT VARIABLES INTEGER ,INTENT(IN) :: nrho !! number of surfaces, (second dimension of LA_in and nu_in modes) TYPE(t_fbase),INTENT(IN) :: fbase_in !< same basis of lambda and nu REAL(wp) ,INTENT(IN) :: LA_in(1:fbase_in%modes,nrho) !< fourier coefficients of thet*=thet+LA(theta,zeta)+iota*nu(theta,zeta) INTEGER ,INTENT(IN) :: tz_dim !< size of the list in thetstar,zetastar REAL(wp) ,INTENT(IN) :: tz_pest(2,tz_dim) !< theta,zeta positions in pest angle (same for all rho) !----------------------------------------------------------------------------------------------------------------------------------- ! OUTPUT VARIABLES REAL(wp) ,INTENT(OUT) :: thetzeta_out(2,tz_dim,nrho) !! theta,zeta position in original angles, for given pest angles !----------------------------------------------------------------------------------------------------------------------------------- ! LOCAL VARIABLES INTEGER :: irho,j REAL(wp) :: theta_star,zeta REAL(wp) :: check(tz_dim),maxerr(nrho) LOGICAL :: docheck !=================================================================================================================================== __PERFON('find_pest_angles') docheck=(testlevel.GT.0) SWRITE(UNIT_StdOut,'(A,2(I8,A))')'Find pest angles via Newton on nrho=',nrho,' times ntheta_zeta= ',tz_dim, " points" CALL ProgressBar(0,nrho)!init DO irho=1,nrho !$OMP PARALLEL DO SCHEDULE(STATIC) DEFAULT(NONE) & !$OMP PRIVATE(j,theta_star,zeta) FIRSTPRIVATE(irho) & !$OMP SHARED(tz_dim,tz_pest,thetzeta_out,fbase_in,LA_in) DO j=1,tz_dim theta_star=tz_pest(1,j); zeta=tz_pest(2,j) thetzeta_out(1,j,irho)=get_pest_newton(theta_star,zeta,fbase_in,LA_in(:,irho)) thetzeta_out(2,j,irho)=zeta END DO! j !$OMP END PARALLEL DO CALL ProgressBar(irho,nrho) END DO !irho IF(docheck)THEN __PERFON("check") DO irho=1,nrho check=fbase_in%evalDOF_xn(tz_dim,thetzeta_out(:,:,irho),0,LA_in(:,irho)) maxerr(irho)=maxval(abs(check+(thetzeta_out(1,:,irho)-tz_pest(1,:)))) END DO IF(ANY(maxerr.GT.1.0e-12))THEN WRITE(UNIT_stdout,*)'CHECK PEST THETA*',maxerr CALL abort(__STAMP__, & "Find_pest_Angles: Error in theta*") END IF __PERFOFF("check") END IF! docheck SWRITE(UNIT_StdOut,'(A)') '...DONE.' __PERFOFF('find_pest_angles') END SUBROUTINE find_pest_angles