38 integer :: n, idir, igrid, ipe_particle, nparticles_local
44 double precision :: w(ixg^t,1:nw)
57 rrd(n,idir) =
rng%unif_01()
61 {^
d&x(^
d,n) = xprobmin^
d + rrd(n+1,^
d) * (xprobmax^
d - xprobmin^
d)\}
78 if(ipe_particle ==
mype)
then
89 nparticles_local = nparticles_local + 1
99 call phys_to_primitive(ixg^
ll,ixg^
ll,w,ps(igrid)%x)
102 w(ixg^t,iw_mom(idir)),ps(igrid)%x,x(:,n),v(idir,n))
108 call advect_update_payload(igrid,ps(igrid)%w,pso(igrid)%w,ps(igrid)%x,x(:,n),v(:,n),q(n),m(n),defpayload,
ndefpayload,0.d0)
111 call usr_update_payload(igrid,ps(igrid)%w,pso(igrid)%w,ps(igrid)%x,x(:,n),v(:,n),q(n),m(n),usrpayload,
nusrpayload,0.d0)
125 integer :: igrid, iigrid, idir
126 double precision,
dimension(ixG^T,1:nw) :: w
128 do iigrid=1,igridstail; igrid=igrids(iigrid);
132 gridvars(igrid)%w(ixg^t,1:
ngridvars) = 0.0d0
133 w(ixg^t,1:nw) = ps(igrid)%w(ixg^t,1:nw)
134 call phys_to_primitive(ixg^
ll,ixg^
ll,w,ps(igrid)%x)
136 gridvars(igrid)%w(ixg^t,
vp(:)) = w(ixg^t,iw_mom(:))
139 gridvars(igrid)%wold(ixg^t,1:
ngridvars) = 0.0d0
140 w(ixg^t,1:nw) = pso(igrid)%w(ixg^t,1:nw)
141 call phys_to_primitive(ixg^
ll,ixg^
ll,w,ps(igrid)%x)
142 gridvars(igrid)%wold(ixg^t,
vp(:)) = w(ixg^t,iw_mom(:))
154 double precision,
intent(in) :: end_time
156 double precision,
dimension(1:ndir) :: v, x
157 double precision :: defpayload(ndefpayload)
158 double precision :: usrpayload(nusrpayload)
159 double precision :: tloc, tlocnew, dt_p, h1
160 double precision,
parameter :: eps=1.0
d-6, hmin=1.0
d-8
161 integer :: ipart, iipart, igrid
162 integer :: nok, nbad, ierror
185 call odeint(x,
ndir,tloc,tlocnew,eps,h1,hmin,nok,nbad,
derivs_advect,
rkqs,ierror)
187 if (ierror /= 0)
then
188 print *,
"odeint returned error code", ierror
189 print *,
"1 means hmin too small, 2 means MAXSTP exceeded"
190 print *,
"Having a problem with particle", iipart
203 call advect_update_payload(igrid,ps(igrid)%w,pso(igrid)%w,ps(igrid)%x,x,v,0.d0,0.d0,defpayload,ndefpayload,tlocnew)
204 particle(ipart)%payload(1:ndefpayload) = defpayload
206 call usr_update_payload(igrid,ps(igrid)%w,pso(igrid)%w,ps(igrid)%x,x,v,0.d0,0.d0,usrpayload,nusrpayload,tlocnew)
215 subroutine advect_update_payload(igrid,w,wold,xgrid,xpart,upart,qpart,mpart,mypayload,mynpayload,particle_time)
217 integer,
intent(in) :: igrid,mynpayload
218 double precision,
intent(in) :: w(ixG^T,1:nw),wold(ixG^T,1:nw)
219 double precision,
intent(in) :: xgrid(ixG^T,1:ndim),xpart(1:ndir),upart(1:ndir),qpart,mpart,particle_time
220 double precision,
intent(out) :: mypayload(mynpayload)
221 double precision :: rho, rho1, rho2, td
226 if (mynpayload > 0 )
then
232 rho = rho1 * (1.0d0 - td) + rho2 * td
241 double precision :: t_s, x(ndir)
242 double precision :: dxdt(ndir)
244 double precision :: v(ndir)
255 double precision,
intent(in) :: end_time
256 double precision :: dt_p
257 integer :: ipart, iipart, nout
258 double precision :: tout, dt_cfl
259 double precision :: v(1:
ndir)
272 dt_cfl = min({
rnode(
rpdx^
d_,partp%igrid)/v(^
d)},bigdouble)
276 if(
phi_ .gt.
ndim) dt_cfl = min(dt_cfl, &
277 sqrt(
rnode(rpdx1_,partp%igrid)/partp%self%x(
r_)) &
280 dt_cfl = min(dt_cfl,0.1d0/v(
phi_))
282 dt_cfl = min(dt_cfl,(partp%self%x(
r_)+smalldouble)/v(
r_))
295 integer,
intent(in) :: igrid, ibeg, iend
296 double precision,
dimension(ndir),
intent(in) :: x
297 double precision,
intent(in) :: tloc
298 double precision,
dimension(iend-ibeg+1),
intent(out) :: var
299 double precision,
dimension(iend-ibeg+1) :: e1, e2
300 integer :: ivar, iloc
301 double precision :: td
314 var(iloc) = e1(iloc) * (1.0d0 - td) + e2(iloc) * td
Module with geometry-related routines (e.g., divergence, curl)
integer, parameter cylindrical
This module contains definitions of global parameters and variables and some generic functions/subrou...
double precision, dimension(:), allocatable w_convert_factor
Conversion factors the primitive variables.
double precision global_time
The global simulation time.
integer, parameter ndim
Number of spatial dimensions for grid variables.
integer icomm
The MPI communicator.
integer mype
The rank of the current MPI task.
integer, dimension(:), allocatable, parameter d
integer ndir
Number of spatial dimensions (components) for vector variables.
integer ixm
the mesh range (within a block with ghost cells)
integer ierrmpi
A global MPI error return code.
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
integer r_
Indices for cylindrical coordinates FOR TESTS, negative value when not used:
double precision, dimension(ndim) dxlevel
This module packages odeint from numerical recipes.
subroutine odeint(ystart, nvar, x1, x2, eps, h1, hmin, nok, nbad, derivs, rkqs, ierror)
subroutine rkqs(y, dydx, n, x, htry, eps, yscal, hdid, hnext, derivs)
Tracer for advected particles moving with fluid flows By Jannis Teunissen, Bart Ripperda,...
subroutine, public advect_init()
pure double precision function advect_get_particle_dt(partp, end_time)
subroutine advect_integrate_particles(end_time)
subroutine advect_fill_gridvars
subroutine derivs_advect(t_s, x, dxdt)
subroutine advect_update_payload(igrid, w, wold, xgrid, xpart, upart, qpart, mpart, mypayload, mynpayload, particle_time)
Payload update.
subroutine get_vec_advect(igrid, x, tloc, var, ibeg, iend)
subroutine, public advect_create_particles()
Module with shared functionality for all the particle movers.
pure subroutine limit_dt_endtime(t_left, dt_p)
integer num_particles
Number of particles.
double precision const_dt_particles
If positive, a constant time step for the particles.
integer ngridvars
Number of variables for grid field.
integer nusrpayload
Number of user-defined payload variables for a particle.
subroutine find_particle_ipe(x, igrid_particle, ipe_particle)
integer npayload
Number of total payload variables for a particle.
integer, dimension(:), allocatable particles_active_on_mype
Array of identity numbers of active particles in current processor.
subroutine push_particle_into_particles_on_mype(ipart)
integer nparticles_active_on_mype
Number of active particles in current processor.
procedure(sub_define_additional_gridvars), pointer particles_define_additional_gridvars
integer nparticles
Identity number and total number of particles.
type(particle_ptr), dimension(:), allocatable particle
integer, dimension(:), allocatable vp
Variable index for fluid velocity.
procedure(sub_integrate), pointer particles_integrate
procedure(sub_fill_gridvars), pointer particles_fill_gridvars
integer ndefpayload
Number of default payload variables for a particle.
integer igrid_working
set the current igrid for the particle integrator:
subroutine interpolate_var(igrid, ixIL, ixOL, gf, x, xloc, gfloc)
Module with all the methods that users can customize in AMRVAC.
procedure(check_particle), pointer usr_check_particle
procedure(create_particles), pointer usr_create_particles
procedure(update_payload), pointer usr_update_payload