12 integer,
private,
parameter :: rho_ = 1
18 character(len=*),
intent(in) :: files(:)
24 open(
unitpar, file=trim(files(n)), status=
"old")
25 read(
unitpar, rotating_frame_list,
end=111)
46 integer,
intent(in) :: ixI^L, ixO^L
47 double precision,
intent(in) :: qdt, x(ixI^S,1:ndim)
48 double precision,
intent(in) :: wCT(ixI^S,1:nw)
49 double precision,
intent(inout) :: w(ixI^S,1:nw)
54 double precision :: rotating_terms(ixI^S),frame_omega(ixI^S)
58 rotating_terms(ixo^s) =
omega_frame**2 * x(ixo^s,
r_) * wct(ixo^s,iw_rho)
60 rotating_terms(ixo^s) = rotating_terms(ixo^s) + 2.d0 *
omega_frame *wct(ixo^s,iw_mom(
phi_))
62 w(ixo^s, iw_mom(
r_)) = w(ixo^s, iw_mom(
r_)) + qdt * rotating_terms(ixo^s)
65 rotating_terms(ixo^s) = - two *
omega_frame * wct(ixo^s,iw_mom(
r_))
66 w(ixo^s, iw_mom(
phi_)) = w(ixo^s, iw_mom(
phi_)) + qdt * rotating_terms(ixo^s)
73 rotating_terms(ixo^s) = frame_omega(ixo^s)**2 * x(ixo^s,
r_) * wct(ixo^s,iw_rho)
75 rotating_terms(ixo^s) = rotating_terms(ixo^s) + &
76 2.d0 * frame_omega(ixo^s) * wct(ixo^s,iw_mom(
phi_))
78 w(ixo^s, iw_mom(
r_)) = w(ixo^s, iw_mom(
r_)) + qdt * rotating_terms(ixo^s)
82 rotating_terms(ixo^s) = frame_omega(ixo^s)**2 * x(ixo^s,
r_) * wct(ixo^s,iw_rho)
84 rotating_terms(ixo^s) = rotating_terms(ixo^s) + &
85 2.d0*wct(ixo^s, iw_mom(
phi_))* frame_omega(ixo^s)
87 w(ixo^s, iw_mom(2)) = w(ixo^s, iw_mom(2)) + qdt * rotating_terms(ixo^s)/ tan(x(ixo^s, 2))
91 rotating_terms(ixo^s) = -2.d0*frame_omega(ixo^s)* wct(ixo^s, iw_mom(
r_))&
92 - 2.d0*wct(ixo^s, iw_mom(2)) * frame_omega(ixo^s)/ tan(x(ixo^s, 2))
93 w(ixo^s, iw_mom(3)) = w(ixo^s, iw_mom(3)) + qdt * rotating_terms(ixo^s)
97 call mpistop(
"Rotating frame not implemented in this geometries")
108 integer,
intent(in) :: ixI^L, ixO^L
109 double precision,
intent(in) :: x(ixI^S,1:ndim)
110 double precision,
intent(out) :: frame_vel(ixI^S)
115 frame_vel(ixo^s) = frame_vel(ixo^s) * dsin(x(ixo^s,2))
126 integer,
intent(in) :: ixI^L, ixO^L
127 double precision,
intent(in) :: x(ixI^S,1:ndim)
128 double precision,
intent(out) :: frame_omega(ixI^S)
133 frame_omega(ixo^s) = frame_omega(ixo^s)* dsin(x(ixo^s,2))
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
Module with geometry-related routines (e.g., divergence, curl)
integer, parameter spherical
integer, parameter cylindrical
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer, parameter unitpar
file handle for IO
logical angmomfix
Enable to strictly conserve the angular momentum (works both in cylindrical and spherical coordinates...
character(len=std_len), dimension(:), allocatable par_files
Which par files are used as input.
integer r_
Indices for cylindrical coordinates FOR TESTS, negative value when not used:
Module for including rotating frame in hydrodynamics simulations The rotation vector is assumed to be...
subroutine rotating_frame_params_read(files)
Read this module's parameters from a file.
subroutine rotating_frame_add_source(qdt, ixIL, ixOL, wCT, w, x)
w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO
double precision omega_frame
Rotation frequency of the frame.
subroutine rotating_frame_omega(x, ixIL, ixOL, frame_omega)
subroutine rotating_frame_init()
Initialize the module.
subroutine rotating_frame_velocity(x, ixIL, ixOL, frame_vel)
Module with all the methods that users can customize in AMRVAC.