-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathsetup_dustydiffusion.f90
126 lines (119 loc) · 3.97 KB
/
setup_dustydiffusion.f90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
!------------------------------------------------------------------------------!
! NDSPMHD: A Smoothed Particle (Magneto)Hydrodynamics code for (astrophysical) !
! fluid dynamics simulations in 1, 2 and 3 spatial dimensions. !
! !
! (c) 2002-2015 Daniel Price !
! !
! http://users.monash.edu.au/~dprice/ndspmhd !
! daniel.price@monash.edu -or- dprice@cantab.net (forwards to current address) !
! !
! NDSPMHD comes with ABSOLUTELY NO WARRANTY. !
! This is free software; and you are welcome to redistribute !
! it under the terms of the GNU General Public License !
! (see LICENSE file for details) and the provision that !
! this notice remains intact. If you modify this file, please !
! note section 2a) of the GPLv2 states that: !
! !
! a) You must cause the modified files to carry prominent notices !
! stating that you changed the files and the date of any change. !
! !
! ChangeLog: !
!------------------------------------------------------------------------------!
!----------------------------------------------------------------
! Set up a uniform density cartesian grid of particles in ND
!----------------------------------------------------------------
subroutine setup
!
!--include relevant global variables
!
use dimen_mhd
use debug
use loguns
use bound
use options
use part
use setup_params
use eos, only:polyk
use uniform_distributions
!
!--define local variables
!
implicit none
integer :: i,iprofile
real :: massp,volume,totmass
real :: denszero,dustfrac0,A,B,r,rmin,rcrit
!
!--allow for tracing flow
!
if (trace) write(iprint,*) ' entering subroutine setup(dustydiffusion)'
!
!--set boundaries
!
ibound = 3 ! boundaries
nbpts = 0 ! use ghosts not fixed
xmin(:) = -0.5 ! set position of boundaries
xmax(:) = 0.5
call set_uniform_cartesian(2,psep,xmin,xmax,adjustbound=.true.)
npart = ntotal
print*,'npart =',npart
if (idust.eq.2) then
itype(1:npart) = itypegas
print*,' setting up dust particles...'
!call set_uniform_cartesian(2,psep,xmin,xmax,adjustbound=.true.)
endif
!
!--determine particle mass
!
denszero = 1.0
volume = product(xmax(:)-xmin(:))
totmass = denszero*volume
massp = totmass/float(ntotal) ! average particle mass
!
!--now assign particle properties
!
dustfrac0 = 0.1 !1.e-3
rcrit = 0.25
B = rcrit**2/dustfrac0
A = dustfrac0/B**(-0.6)
print*,' A = ',A,' B = ',B
iprofile = 2
polyk = 1.
rmin = psep
do i=1,ntotal
vel(:,i) = 0.
dens(i) = denszero
pmass(i) = massp
uu(i) = 1.5*polyk ! isothermal
Bfield(:,i) = 0.
if (onef_dust) then
r = sqrt(dot_product(x(:,i),x(:,i)))
select case(iprofile)
case(3)
dustfrac(:,i) = 0.5*exp(-(r/0.1)**2)
case(2)
dustfrac(:,i) = max(dustfrac0*(1. - (r/rcrit)**2),1.e-4)
case(1)
dustfrac(:,i) = r**2/A
case default
if (r < 0.5) then
dustfrac(:,i) = min(sqrt(A**2/r + dustfrac0**2),0.9)
else
dustfrac(:,i) = dustfrac0
endif
end select
else
dustfrac(:,i) = 0.
endif
enddo
!
!--allow for tracing flow
!
if (trace) write(iprint,*) ' exiting subroutine setup'
return
end
!
! use this routine to modify the dump upon code restart
!
subroutine modify_dump()
implicit none
end subroutine modify_dump