Skip to content

Commit

Permalink
Merge pull request #206 from PrincetonUniversity/bugfix/BEAMS3D_minor
Browse files Browse the repository at this point in the history
Bugfix/beams3 d minor
  • Loading branch information
lazersos authored Oct 26, 2023
2 parents d374529 + 631f90e commit af17e3a
Show file tree
Hide file tree
Showing 7 changed files with 317 additions and 255 deletions.
1 change: 1 addition & 0 deletions BEAMS3D/BEAMS3D.dep
Original file line number Diff line number Diff line change
Expand Up @@ -366,6 +366,7 @@ fidasim_input_mod.o: \
../../LIBSTELL/$(LOCTYPE)/wall_mod.o \
../../LIBSTELL/$(LOCTYPE)/mpi_params.o \
../../LIBSTELL/$(LOCTYPE)/mpi_inc.o \
beams3d_physics_mod.o \
beams3d_runtime.o \
beams3d_lines.o \
beams3d_write_parhdf5.o \
Expand Down
4 changes: 2 additions & 2 deletions BEAMS3D/Sources/beams3d_diagnostics.f90
Original file line number Diff line number Diff line change
Expand Up @@ -107,8 +107,8 @@ SUBROUTINE beams3d_diagnostics
CALL MPI_REDUCE(shine_port, shine_port, nbeams, MPI_DOUBLE_PRECISION, MPI_SUM, master, MPI_COMM_BEAMS, ierr_mpi)
CALL MPI_REDUCE(norbit, norbit, nbeams, MPI_REAL, MPI_SUM, master, MPI_COMM_BEAMS, ierr_mpi)
CALL MPI_REDUCE(nlost, nlost, nbeams, MPI_REAL, MPI_SUM, master, MPI_COMM_BEAMS, ierr_mpi)
CALL MPI_REDUCE(nlost, tlow, nbeams, MPI_REAL, MPI_MIN, master, MPI_COMM_BEAMS, ierr_mpi)
CALL MPI_REDUCE(nlost, thigh, nbeams, MPI_REAL, MPI_MAX, master, MPI_COMM_BEAMS, ierr_mpi)
CALL MPI_REDUCE(tlow, tlow, nbeams, MPI_REAL, MPI_MIN, master, MPI_COMM_BEAMS, ierr_mpi)
CALL MPI_REDUCE(thigh, thigh, nbeams, MPI_REAL, MPI_MAX, master, MPI_COMM_BEAMS, ierr_mpi)
END IF
#endif

Expand Down
5 changes: 4 additions & 1 deletion BEAMS3D/Sources/beams3d_follow.f90
Original file line number Diff line number Diff line change
Expand Up @@ -302,11 +302,14 @@ SUBROUTINE beams3d_follow
DEALLOCATE(itemp)
IF (ALLOCATED(mnum)) DEALLOCATE(mnum)
IF (ALLOCATED(moffsets)) DEALLOCATE(moffsets)
IF (ALLOCATED(t_last)) DEALLOCATE(t_last)
CALL MPI_BARRIER(MPI_COMM_BEAMS, ierr_mpi)
IF (ierr_mpi /= 0) CALL handle_err(MPI_BARRIER_ERR, 'beams3d_follow', ierr_mpi)
#endif

! Adjust T_END back to values of T_last
t_end(mystart:myend) = t_last(mystart:myend)
IF (ALLOCATED(t_last)) DEALLOCATE(t_last)

RETURN
!-----------------------------------------------------------------------
! End Subroutine
Expand Down
2 changes: 1 addition & 1 deletion BEAMS3D/Sources/beams3d_init.f90
Original file line number Diff line number Diff line change
Expand Up @@ -423,7 +423,7 @@ SUBROUTINE beams3d_init
! Load vessel if not done already vessel
IF (lvessel .and. (.not. lwall_loaded)) THEN
IF (lverb) THEN
WRITE(6,'(A)') 'Loading Vessel!'
WRITE(6,'(A)') '----- Loading wall data -----'
END IF
CALL wall_load_txt(TRIM(vessel_string),ier,.false.,MPI_COMM_BEAMS)
IF (lverb) THEN
Expand Down
13 changes: 8 additions & 5 deletions BEAMS3D/Sources/beams3d_interface_mod.f90
Original file line number Diff line number Diff line change
Expand Up @@ -78,12 +78,15 @@ SUBROUTINE beams3d_cleanup
CALL beams3d_free(MPI_COMM_SHARMEM)
IF (lvessel) CALL wall_free(ier,MPI_COMM_BEAMS)
#if defined(MPI_OPT)
CALL MPI_BARRIER(MPI_COMM_BEAMS, ierr_mpi)
ierr_mpi = 0; CALL MPI_BARRIER(MPI_COMM_BEAMS, ierr_mpi)
IF (ierr_mpi /= 0) CALL handle_err(MPI_BARRIER_ERR, 'beams3d_main', ierr_mpi)
ierr_mpi=0
CALL MPI_INFO_FREE(mpi_info_beams3d, ierr_mpi)
CALL MPI_FINALIZE(ierr_mpi)
IF (ierr_mpi /= MPI_SUCCESS) CALL handle_err(MPI_FINE_ERR, 'beams3d_main', ierr_mpi)
ierr_mpi = 0; CALL MPI_INFO_FREE(mpi_info_beams3d, ierr_mpi)
ierr_mpi = 0; CALL MPI_BARRIER(MPI_COMM_BEAMS, ierr_mpi)
IF (ierr_mpi /= 0) CALL handle_err(MPI_BARRIER_ERR, 'beams3d_main', ierr_mpi)
ierr_mpi = 0; CALL MPI_COMM_FREE(MPI_COMM_SHARMEM, ierr_mpi)
ierr_mpi = 0; CALL MPI_COMM_FREE(MPI_COMM_BEAMS, ierr_mpi)
ierr_mpi = 0; CALL MPI_FINALIZE(ierr_mpi)
!IF (ierr_mpi /= MPI_SUCCESS) CALL handle_err(MPI_FINE_ERR, 'beams3d_main', ierr_mpi)
#endif
IF (lverb) WRITE(6, '(A)') '----- BEAMS3D DONE -----'
RETURN
Expand Down
72 changes: 72 additions & 0 deletions BEAMS3D/Sources/beams3d_physics_mod.f90
Original file line number Diff line number Diff line change
Expand Up @@ -1811,6 +1811,78 @@ SUBROUTINE beams3d_SFLX(q,S)

END SUBROUTINE beams3d_SFLX

!-----------------------------------------------------------------
! Function: beams3d_SUFLX
! Authors: S. Lazerson (samuel.lazerson@ipp.mpg.de)
! Date: 09/28/2023
! Description: Returns normalized toroidal flux and
! poloidal angle at a given point in space.
!-----------------------------------------------------------------
SUBROUTINE beams3d_SUFLX(q,S,U)
!--------------------------------------------------------------
! Input Parameters
! q (q(1),q(2),q(3)) = (R,phi,Z)
! S Backbround grid normalized flux
! U Backbround grid poloidal angle
!--------------------------------------------------------------
IMPLICIT NONE
DOUBLE PRECISION, INTENT(inout) :: q(3)
DOUBLE PRECISION, INTENT(out) :: S
DOUBLE PRECISION, INTENT(out) :: U

!--------------------------------------------------------------
! Local Variables
! r_temp Helpers (r,phi,z)
! i,j,k Spline Grid indicies
! xparam Spline subgrid factor [0,1] (yparam,zparam)
! ict Spline output control
! fval Spline output array
!--------------------------------------------------------------
DOUBLE PRECISION :: r_temp, z_temp, phi_temp
! For splines
INTEGER :: i,j,k
REAL*8 :: xparam, yparam, zparam
INTEGER, parameter :: ict(8)=(/1,0,0,0,0,0,0,0/)
REAL*8 :: fval(1)

!--------------------------------------------------------------
! Begin Subroutine
!--------------------------------------------------------------

! Setup position in a vll arrays
r_temp = q(1)
phi_temp = MODULO(q(2), phimax)
IF (phi_temp < 0) phi_temp = phi_temp + phimax
z_temp = q(3)

! Initialize values
S = 2; U = 0

! Check that we're inside the domain then proceed
IF ((r_temp >= rmin-eps1) .and. (r_temp <= rmax+eps1) .and. &
(phi_temp >= phimin-eps2) .and. (phi_temp <= phimax+eps2) .and. &
(z_temp >= zmin-eps3) .and. (z_temp <= zmax+eps3)) THEN
i = MIN(MAX(COUNT(raxis < r_temp),1),nr-1)
j = MIN(MAX(COUNT(phiaxis < phi_temp),1),nphi-1)
k = MIN(MAX(COUNT(zaxis < z_temp),1),nz-1)
xparam = (r_temp - raxis(i)) * hri(i)
yparam = (phi_temp - phiaxis(j)) * hpi(j)
zparam = (z_temp - zaxis(k)) * hzi(k)
! Evaluate the Splines
CALL R8HERM3FCN(ict,1,1,fval,i,j,k,xparam,yparam,zparam,&
hr(i),hri(i),hp(j),hpi(j),hz(k),hzi(k),&
S4D(1,1,1,1),nr,nphi,nz)
S = max(fval(1),zero)
CALL R8HERM3FCN(ict,1,1,fval,i,j,k,xparam,yparam,zparam,&
hr(i),hri(i),hp(j),hpi(j),hz(k),hzi(k),&
U4D(1,1,1,1),nr,nphi,nz)
U = fval(1)
END IF

RETURN

END SUBROUTINE beams3d_SUFLX

!-----------------------------------------------------------------
! Function: beams3d_suv2rzp
! Authors: S. Lazerson (samuel.lazerson@ipp.mpg.de)
Expand Down
Loading

0 comments on commit af17e3a

Please sign in to comment.