Skip to content
This repository was archived by the owner on Aug 28, 2024. It is now read-only.

Commit

Permalink
Created basic tidal acceleration model. Currently there is no feedbac…
Browse files Browse the repository at this point in the history
…k between the tide and oblateness acceleration models, so this should be done in the near future once the spin integration is done
  • Loading branch information
daminton committed Jul 13, 2021
1 parent 3204af0 commit 1181c9a
Show file tree
Hide file tree
Showing 10 changed files with 97 additions and 937 deletions.
6 changes: 5 additions & 1 deletion src/helio/helio_getacch.f90
Original file line number Diff line number Diff line change
Expand Up @@ -22,9 +22,13 @@ module subroutine helio_getacch_pl(self, system, param, t, lbeg)
cb%aoblbeg = cb%aobl
call pl%accel_obl(system)
cb%aoblend = cb%aobl
if (param%ltides) then
cb%atidebeg = cb%atide
call pl%accel_tides(system)
cb%atideend = cb%atide
end if
end if
if (param%lextra_force) call pl%accel_user(system, param, t)
if (param%ltides) call pl%accel_tides(system)
!if (param%lgr) call pl%gr_accel(param)
end associate

Expand Down
7 changes: 4 additions & 3 deletions src/modules/rmvs_classes.f90
Original file line number Diff line number Diff line change
Expand Up @@ -31,9 +31,10 @@ module rmvs_classes
end type rmvs_nbody_system

type, private :: rmvs_interp
real(DP), dimension(:, :), allocatable :: x !! interpolated heliocentric planet position for outer encounter
real(DP), dimension(:, :), allocatable :: v !! interpolated heliocentric planet velocity for outer encounter
real(DP), dimension(:, :), allocatable :: aobl !! Encountering planet's oblateness acceleration value
real(DP), dimension(:, :), allocatable :: x !! interpolated heliocentric planet position for outer encounter
real(DP), dimension(:, :), allocatable :: v !! interpolated heliocentric planet velocity for outer encounter
real(DP), dimension(:, :), allocatable :: aobl !! Encountering planet's oblateness acceleration value
real(DP), dimension(:, :), allocatable :: atide !! Encountering planet's tidal acceleration value
end type rmvs_interp

!********************************************************************************************************************************
Expand Down
54 changes: 25 additions & 29 deletions src/modules/swiftest_classes.f90
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@ module swiftest_classes
public :: tides_getacch_pl, tides_step_spin_system
public :: user_getacch_body
public :: util_coord_b2h_pl, util_coord_b2h_tp, util_coord_h2b_pl, util_coord_h2b_tp, util_exit, util_fill_body, util_fill_pl, util_fill_tp, &
util_index, util_peri_tp, util_reverse_status, util_set_beg_end_cb, util_set_beg_end_pl, util_set_ir3h, util_set_msys, util_set_mu_pl, &
util_index, util_peri_tp, util_reverse_status, util_set_beg_end_pl, util_set_ir3h, util_set_msys, util_set_mu_pl, &
util_set_mu_tp, util_set_rhill, util_sort, util_spill_body, util_spill_pl, util_spill_tp, util_valid, util_version

!********************************************************************************************************************************
Expand Down Expand Up @@ -104,32 +104,34 @@ module swiftest_classes
!> A concrete lass for the central body in a Swiftest simulation
type, abstract, public, extends(swiftest_base) :: swiftest_cb
character(len=STRMAX) :: name !! Non-unique name
integer(I4B) :: id = 0 !! External identifier (unique)
real(DP) :: mass = 0.0_DP !! Central body mass (units MU)
real(DP) :: Gmass = 0.0_DP !! Central mass gravitational term G * mass (units GU * MU)
real(DP) :: radius = 0.0_DP !! Central body radius (units DU)
real(DP) :: density = 1.0_DP !! Central body mass density - calculated internally (units MU / DU**3)
real(DP) :: j2rp2 = 0.0_DP !! J2*R^2 term for central body
real(DP) :: j4rp4 = 0.0_DP !! J4*R^2 term for central body
real(DP), dimension(NDIM) :: aobl = 0.0_DP !! Barycentric acceleration due to central body oblatenes
real(DP), dimension(NDIM) :: aoblbeg = 0.0_DP !! Barycentric acceleration due to central body oblatenes at beginning of step
real(DP), dimension(NDIM) :: aoblend = 0.0_DP !! Barycentric acceleration due to central body oblatenes at end of step
real(DP), dimension(NDIM) :: xb = 0.0_DP !! Barycentric position (units DU)
real(DP), dimension(NDIM) :: vb = 0.0_DP !! Barycentric velocity (units DU / TU)
real(DP), dimension(NDIM) :: agr = 0.0_DP !! Acceleration due to post-Newtonian correction
real(DP), dimension(NDIM) :: Ip = 0.0_DP !! Unitless principal moments of inertia (I1, I2, I3) / (MR**2). Principal axis rotation assumed.
real(DP), dimension(NDIM) :: rot = 0.0_DP !! Body rotation vector in inertial coordinate frame (units rad / TU)
real(DP) :: k2 = 0.0_DP !! Tidal Love number
real(DP) :: Q = 0.0_DP !! Tidal quality factor
real(DP) :: tlag = 0.0_DP !! Tidal phase lag angle
real(DP), dimension(NDIM) :: L0 = 0.0_DP !! Initial angular momentum of the central body
real(DP), dimension(NDIM) :: dL = 0.0_DP !! Change in angular momentum of the central body
integer(I4B) :: id = 0 !! External identifier (unique)
real(DP) :: mass = 0.0_DP !! Central body mass (units MU)
real(DP) :: Gmass = 0.0_DP !! Central mass gravitational term G * mass (units GU * MU)
real(DP) :: radius = 0.0_DP !! Central body radius (units DU)
real(DP) :: density = 1.0_DP !! Central body mass density - calculated internally (units MU / DU**3)
real(DP) :: j2rp2 = 0.0_DP !! J2*R^2 term for central body
real(DP) :: j4rp4 = 0.0_DP !! J4*R^2 term for central body
real(DP), dimension(NDIM) :: aobl = 0.0_DP !! Barycentric acceleration due to central body oblatenes
real(DP), dimension(NDIM) :: atide = 0.0_DP !! Barycentric acceleration due to central body oblatenes
real(DP), dimension(NDIM) :: aoblbeg = 0.0_DP !! Barycentric acceleration due to central body oblatenes at beginning of step
real(DP), dimension(NDIM) :: aoblend = 0.0_DP !! Barycentric acceleration due to central body oblatenes at end of step
real(DP), dimension(NDIM) :: atidebeg = 0.0_DP !! Barycentric acceleration due to central body oblatenes at beginning of step
real(DP), dimension(NDIM) :: atideend = 0.0_DP !! Barycentric acceleration due to central body oblatenes at end of step
real(DP), dimension(NDIM) :: xb = 0.0_DP !! Barycentric position (units DU)
real(DP), dimension(NDIM) :: vb = 0.0_DP !! Barycentric velocity (units DU / TU)
real(DP), dimension(NDIM) :: agr = 0.0_DP !! Acceleration due to post-Newtonian correction
real(DP), dimension(NDIM) :: Ip = 0.0_DP !! Unitless principal moments of inertia (I1, I2, I3) / (MR**2). Principal axis rotation assumed.
real(DP), dimension(NDIM) :: rot = 0.0_DP !! Body rotation vector in inertial coordinate frame (units rad / TU)
real(DP) :: k2 = 0.0_DP !! Tidal Love number
real(DP) :: Q = 0.0_DP !! Tidal quality factor
real(DP) :: tlag = 0.0_DP !! Tidal phase lag angle
real(DP), dimension(NDIM) :: L0 = 0.0_DP !! Initial angular momentum of the central body
real(DP), dimension(NDIM) :: dL = 0.0_DP !! Change in angular momentum of the central body
contains
private
procedure, public :: initialize => io_read_cb_in !! I/O routine for reading in central body data
procedure, public :: write_frame => io_write_frame_cb !! I/O routine for writing out a single frame of time-series data for the central body
procedure, public :: read_frame => io_read_frame_cb !! I/O routine for reading out a single frame of time-series data for the central body
procedure, public :: set_beg_end => util_set_beg_end_cb !! Sets the beginning and ending oblateness acceleration term
end type swiftest_cb

!********************************************************************************************************************************
Expand All @@ -150,6 +152,7 @@ module swiftest_classes
real(DP), dimension(:,:), allocatable :: vb !! Barycentric velocity
real(DP), dimension(:,:), allocatable :: ah !! Total heliocentric acceleration
real(DP), dimension(:,:), allocatable :: aobl !! Barycentric accelerations of bodies due to central body oblatenes
real(DP), dimension(:,:), allocatable :: atide !! Tanngential component of acceleration of bodies due to tides
real(DP), dimension(:,:), allocatable :: agr !! Acceleration due to post-Newtonian correction
real(DP), dimension(:), allocatable :: ir3h !! Inverse heliocentric radius term (1/rh**3)
real(DP), dimension(:), allocatable :: a !! Semimajor axis (pericentric distance for a parabolic orbit)
Expand Down Expand Up @@ -770,13 +773,6 @@ module subroutine util_reverse_status(self)
class(swiftest_body), intent(inout) :: self !! Swiftest body object
end subroutine util_reverse_status

module subroutine util_set_beg_end_cb(self, aoblbeg, aoblend)
implicit none
class(swiftest_cb), intent(inout) :: self !! Swiftest central body object
real(DP), dimension(:), intent(in), optional :: aoblbeg !! Oblateness acceleration term at beginning of step
real(DP), dimension(:), intent(in), optional :: aoblend !! Oblateness acceleration term at end of step
end subroutine util_set_beg_end_cb

module subroutine util_set_beg_end_pl(self, xbeg, xend, vbeg)
implicit none
class(swiftest_pl), intent(inout) :: self !! Swiftest massive body object
Expand Down
Loading

0 comments on commit 1181c9a

Please sign in to comment.