! $Id: particles_density.f90 20849 2013-08-06 18:45:43Z anders@astro.lu.se $ ! ! This module takes care of everything related to the density represented by ! each (super)particle. ! !** AUTOMATIC CPARAM.INC GENERATION **************************** ! ! Declare (for generation of cparam.inc) the number of f array ! variables and auxiliary variables added by this module ! ! MPVAR CONTRIBUTION 1 ! CPARAM logical, parameter :: lparticles_density=.true. ! !*************************************************************** module Particles_density ! use Cdata use General, only: keep_compiler_quiet use Messages use Particles_cdata use Particles_sub ! implicit none ! include 'particles_density.h' ! real :: rhop_swarm0=1.0, rhop_swarm1=1.0, rhop_swarm2, rhop_swarm3 real :: gravr_swarm0=1.0, gravr_swarm1=1.0 real, pointer :: rhs_poisson_const character (len=labellen), dimension(ninit) :: initrhopswarm='nothing' ! namelist /particles_dens_init_pars/ & initrhopswarm, rhop_swarm0, rhop_swarm1, rhop_swarm2, rhop_swarm3, & gravr_swarm0, gravr_swarm1, eps_dtog ! contains !*********************************************************************** subroutine register_particles_density() ! ! Set up indices for access to the fp and dfp arrays. ! ! 22-nov-10/anders+michiel: adapted ! if (lroot) call svn_id( & "$Id: particles_density.f90 20849 2013-08-06 18:45:43Z anders@astro.lu.se $") ! ! Index for particle density. ! call append_npvar('irhopswarm',irhopswarm) ! endsubroutine register_particles_density !*********************************************************************** subroutine initialize_particles_density(f) ! ! Perform any post-parameter-read initialization i.e. calculate derived ! parameters. ! ! 22-nov-10/anders+michiel: adapted ! use SharedVariables, only: get_shared_variable ! real, dimension (mx,my,mz,mfarray) :: f ! if (lselfgravity) then call get_shared_variable('rhs_poisson_const',rhs_poisson_const) endif ! if (.not.(lcartesian_coords.and.(all(lequidist)))) call fatal_error( & 'initialize_particles_density', 'particles_density only implemented '// & 'for Cartesian equidistant grids.') ! call keep_compiler_quiet(f) ! endsubroutine initialize_particles_density !*********************************************************************** subroutine init_particles_density(f,fp) ! ! Initial particle density. ! ! 22-nov-10/anders+michiel: adapted ! real, dimension (mx,my,mz,mfarray) :: f real, dimension (mpar_loc,mparray) :: fp ! real :: rhom integer :: j, k ! do j=1,ninit ! select case (initrhopswarm(j)) ! case ('nothing') if (lroot.and.j==1) print*, 'init_particles_density: nothing' ! case ('constant') if (lroot) then print*, 'init_particles_density: constant particle density' print*, 'init_particles_density: rhop_swarm0=', rhop_swarm0 endif fp(1:npar_loc,irhopswarm)=rhop_swarm0 ! case ('constant-1') if (lroot) then print*, 'init_particles_density: set particle 1 density' print*, 'init_particles_density: rhop_swarm1=', rhop_swarm1 endif do k=1,npar_loc if (ipar(k)==1) fp(k,irhopswarm)=rhop_swarm1 enddo ! case ('constant-2') if (lroot) then print*, 'init_particles_density: set particle 2 density' print*, 'init_particles_density: rhop_swarm2=', rhop_swarm2 endif do k=1,npar_loc if (ipar(k)==2) fp(k,irhopswarm)=rhop_swarm2 enddo ! case ('constant-3') if (lroot) then print*, 'init_particles_density: set particle 3 density' print*, 'init_particles_density: rhop_swarm3=', rhop_swarm3 endif do k=1,npar_loc if (ipar(k)==3) fp(k,irhopswarm)=rhop_swarm3 enddo ! case ('constant-rhop') fp(1:npar_loc,irhopswarm)=rhop_swarm0/(float(npar)/nwgrid) ! case ('constant-grav') if (lroot) then print*, 'init_particles_density: constant particle gravity' print*, 'init_particles_density: gravr_swarm=', gravr_swarm0 endif if (.not. lselfgravity) then if (lroot) print*, 'init_particles_density: need selfgravity '// & 'module for this initial condition' call fatal_error('init_particles_density','') endif fp(1:npar_loc,irhopswarm)=gravr_swarm0/ & (rhs_poisson_const/(4*pi)*dx**3) ! case ('constant-grav-1') if (lroot) then print*, 'init_particles_density: set particle 1 gravity' print*, 'init_particles_density: gravr_swarm1=', gravr_swarm1 endif if (.not. lselfgravity) then if (lroot) print*, 'init_particles_density: need selfgravity '// & 'module for this initial condition' call fatal_error('init_particles_density','') endif do k=1,npar_loc if (ipar(k)==1) & fp(k,irhopswarm)=gravr_swarm1/(rhs_poisson_const/(4*pi)*dx**3) enddo ! case ('from-particles-module','particles-to-gas-ratio') fp(1:npar_loc,irhopswarm)=rhop_swarm ! case ('thin-disk') rhom=sqrt(2*pi)*1.0*1.0/Lz ! rhom = Sigma/Lz, Sigma=sqrt(2*pi)*H*rho1 fp(1:npar_loc,irhopswarm)=eps_dtog*rhom/(real(npar)/nwgrid) ! case default if (lroot) print*, 'init_particles_density: '// & 'No such such value for initrhopswarm: ', trim(initrhopswarm(j)) call fatal_error('init_particles_density','') endselect ! enddo ! call keep_compiler_quiet(f) ! endsubroutine init_particles_density !*********************************************************************** subroutine pencil_criteria_par_density() ! ! All pencils that the Particles_density module depends on are specified ! here. ! ! 22-nov-10/anders+michiel: adapted ! endsubroutine pencil_criteria_par_density !*********************************************************************** subroutine drhopswarm_dt_pencil(f,df,fp,dfp,p,ineargrid) ! ! Evolution of particle density. ! ! 22-nov-10/anders+michiel: adapted ! real, dimension (mx,my,mz,mfarray) :: f real, dimension (mx,my,mz,mvar) :: df real, dimension (mpar_loc,mparray) :: fp real, dimension (mpar_loc,mpvar) :: dfp type (pencil_case) :: p integer, dimension (mpar_loc,3) :: ineargrid ! logical :: lheader, lfirstcall=.true. ! intent (in) :: f, fp intent (out) :: dfp ! ! Print out header information in first time step. ! lheader=lfirstcall .and. lroot ! ! Identify module and boundary conditions. ! if (lheader) print*, 'drhopswarm_dt_pencil: Calculate drhopswarm_dt' ! lfirstcall=.false. ! call keep_compiler_quiet(f) call keep_compiler_quiet(df) call keep_compiler_quiet(fp) call keep_compiler_quiet(dfp) call keep_compiler_quiet(p) call keep_compiler_quiet(ineargrid) ! endsubroutine drhopswarm_dt_pencil !*********************************************************************** subroutine drhopswarm_dt(f,df,fp,dfp,ineargrid) ! ! Evolution of internal particle number. ! ! 22-nov-10/anders+michiel: adapted ! real, dimension (mx,my,mz,mfarray) :: f real, dimension (mx,my,mz,mvar) :: df real, dimension (mpar_loc,mparray) :: fp real, dimension (mpar_loc,mpvar) :: dfp integer, dimension (mpar_loc,3) :: ineargrid ! call keep_compiler_quiet(f,df) call keep_compiler_quiet(fp) call keep_compiler_quiet(dfp) call keep_compiler_quiet(ineargrid) ! endsubroutine drhopswarm_dt !*********************************************************************** subroutine read_particles_dens_init_pars(iostat) ! use File_io, only: parallel_unit ! integer, intent(out) :: iostat ! read(parallel_unit, NML=particles_dens_init_pars, IOSTAT=iostat) ! endsubroutine read_particles_dens_init_pars !*********************************************************************** subroutine write_particles_dens_init_pars(unit) ! integer, intent(in) :: unit ! write(unit, NML=particles_dens_init_pars) ! endsubroutine write_particles_dens_init_pars !*********************************************************************** subroutine rprint_particles_density(lreset,lwrite) ! ! Read and register print parameters relevant for particle density. ! ! 22-nov-10/anders+michiel: adapted ! logical :: lreset logical, optional :: lwrite ! integer :: iname ! call keep_compiler_quiet(lreset) call keep_compiler_quiet(lwrite) ! endsubroutine rprint_particles_density !*********************************************************************** endmodule Particles_density