28 integer,
intent(in) :: ixI^L,ixO^L,idim
29 double precision,
dimension(ixI^S,1:nw),
intent(in) :: wRC,wLC
30 double precision,
dimension(ixI^S,1:nwflux),
intent(in) :: fLC, fRC
31 integer,
dimension(ixI^S),
intent(inout) :: patchf
33 integer :: ixOO^D,TxOO^L
38 {
do ixoo^d= ixo^lim^d\}
43 if(abs(patchf(ixoo^d)) == 1 .or. abs(patchf(ixoo^d)) == 4)
Then
44 if(any(frc(ixoo^d,1:nwflux)*flc(ixoo^d,1:nwflux)<-smalldouble))
Then
45 where(abs(patchf(txoo^s))==1)
54 subroutine mhd_get_lcd(wLC,wRC,fLC,fRC,cmin,cmax,idim,ixI^L,ixO^L, &
55 whll,Fhll,lambdaCD,patchf)
61 integer,
intent(in) :: ixI^L,ixO^L,idim
62 double precision,
dimension(ixI^S,1:nw),
intent(in) :: wLC,wRC
63 double precision,
dimension(ixI^S,1:nwflux),
intent(in) :: fLC,fRC
64 double precision,
dimension(ixI^S),
intent(in) :: cmax,cmin
65 integer ,
dimension(ixI^S),
intent(inout) :: patchf
66 double precision,
dimension(ixI^S,1:nwflux),
intent(out) :: Fhll,whll
67 double precision,
dimension(ixI^S),
intent(out) :: lambdaCD
69 logical ,
dimension(ixO^S) :: Cond_patchf
70 double precision :: Epsilon
79 cond_patchf(ixo^s)=(abs(patchf(ixo^s))==1)
83 where(cond_patchf(ixo^s))
85 fhll(ixo^s,iw)= (cmax(ixo^s)*flc(ixo^s,iw)-cmin(ixo^s)*frc(ixo^s,iw) &
86 + cmin(ixo^s)*cmax(ixo^s)*(wrc(ixo^s,iw)-wlc(ixo^s,iw)))&
87 /(cmax(ixo^s)-cmin(ixo^s))
89 whll(ixo^s,iw) = (cmax(ixo^s)*wrc(ixo^s,iw)-cmin(ixo^s)*wlc(ixo^s,iw)&
90 +flc(ixo^s,iw)-frc(ixo^s,iw))/(cmax(ixo^s)-cmin(ixo^s))
95 where(cond_patchf(ixo^s))
96 lambdacd(ixo^s)=whll(ixo^s,
mom(idim))/whll(ixo^s,
rho_)
99 {
do ix^db=ixomin^db,ixomax^db\}
100 if(cond_patchf(ix^d))
then
103 if(cmin(ix^d)<zero.and.lambdacd(ix^d)>zero&
104 .and.lambdacd(ix^d)<cmax(ix^d))
then
106 else if(cmax(ix^d)>zero.and.lambdacd(ix^d)<zero&
107 .and.lambdacd(ix^d)>cmin(ix^d))
then
109 else if(lambdacd(ix^d)>=cmax(ix^d).or.lambdacd(ix^d) <= cmin(ix^d))
then
110 lambdacd(ix^d) = zero
117 where(patchf(ixo^s)== 3)
118 cond_patchf(ixo^s)=.false.
122 if(any(lambdacd(ixo^s)==zero.and.cond_patchf(ixo^s)))
then
125 {
do ix^db=ixomin^db,ixomax^db\}
126 if(lambdacd(ix^d)==zero.and.cond_patchf(ix^d))
then
127 if(-cmin(ix^d)>=cmax(ix^d))
then
139 where(cond_patchf(ixo^s).and. &
140 dabs(lambdacd(ixo^s))/max(cmax(ixo^s),epsilon)< epsilon .and. &
141 dabs(lambdacd(ixo^s))/max(dabs(cmin(ixo^s)),epsilon)< epsilon)
142 lambdacd(ixo^s) = zero
148 subroutine mhd_get_wcd(wLC,wRC,whll,fRC,fLC,Fhll,patchf,lambdaCD,cmin,cmax,&
157 integer,
intent(in) :: ixI^L,ixO^L,idim
158 double precision,
dimension(ixI^S,1:nw),
intent(in) :: wRC,wLC
159 double precision,
dimension(ixI^S,1:nwflux),
intent(in):: whll, Fhll
160 double precision,
dimension(ixI^S),
intent(in) :: lambdaCD
161 double precision,
dimension(ixI^S),
intent(in) :: cmax,cmin
162 double precision,
dimension(ixI^S,1:nwflux),
intent(in):: fRC,fLC
163 double precision,
dimension(ixI^S,1:nwflux),
intent(out):: f
164 double precision,
dimension(ixI^S,1:nw) :: wCD,wSub
165 double precision,
dimension(ixI^S,1:nwflux) :: fSub
166 double precision,
dimension(ixI^S) :: vSub,cspeed,pCD,VdotBCD
167 integer ,
dimension(ixI^S),
intent(in) :: patchf
169 integer :: n, iw, idir,ix^D
172 {
do ix^db=ixomin^db,ixomax^db\}
173 if(patchf(ix^d)==1)
then
174 cspeed(ix^d)=cmax(ix^d)
175 vsub(ix^d)=wrc(ix^d,mom(idim))/wrc(ix^d,rho_)
176 wsub(ix^d,:)=wrc(ix^d,:)
177 fsub(ix^d,:)=frc(ix^d,:)
178 else if(patchf(ix^d)==-1)
then
179 cspeed(ix^d)=cmin(ix^d)
180 vsub(ix^d)=wlc(ix^d,mom(idim))/wlc(ix^d,rho_)
181 wsub(ix^d,:)=wlc(ix^d,:)
182 fsub(ix^d,:)=flc(ix^d,:)
186 {
do ix^db=ixomin^db,ixomax^db\}
187 if(abs(patchf(ix^d))==1)
then
188 wcd(ix^d,rho_) = wsub(ix^d,rho_)&
189 *(cspeed(ix^d)-vsub(ix^d))/(cspeed(ix^d)-lambdacd(ix^d))
192 wcd(ix^d,iw) = wsub(ix^d,iw)*(cspeed(ix^d)-vsub(ix^d))&
193 /(cspeed(ix^d)-lambdacd(ix^d))
198 wcd(ix^d,mag(idir)) = whll(ix^d,mag(idir))
204 wcd(ix^d,mom(iw))=(cspeed(ix^d)*wsub(ix^d,mom(iw))-fsub(ix^d,mom(iw))&
205 -wcd(ix^d,mag(idim))*wcd(ix^d,mag(iw)))/&
206 (cspeed(ix^d)-lambdacd(ix^d))
209 wcd(ix^d,mom(iw)) = wcd(ix^d,rho_) * lambdacd(ix^d)
213 vdotbcd(ix^d) = sum(whll(ix^d,mom(:))*whll(ix^d,mag(:)))/whll(ix^d,rho_)
215 pcd(ix^d) = wsub(ix^d,rho_)*(cspeed(ix^d)-vsub(ix^d))&
216 *(lambdacd(ix^d)-vsub(ix^d))&
217 +fsub(ix^d,mom(idim))-wsub(ix^d,mom(idim))*vsub(ix^d)&
218 + wcd(ix^d,mag(idim))**2
220 wcd(ix^d,e_) = (cspeed(ix^d)*wsub(ix^d,e_) &
221 -fsub(ix^d,e_)+lambdacd(ix^d)*pcd(ix^d)&
222 -vdotbcd(ix^d)*wcd(ix^d,mag(idim)))&
223 /(cspeed(ix^d)-lambdacd(ix^d))
229 if(iw == mag(idim))
then
231 else if(mhd_glm .and. iw == psi_)
then
234 where(abs(patchf(ixo^s))==1)
236 f(ixo^s,iw)=fsub(ixo^s,iw)+cspeed(ixo^s)*(wcd(ixo^s,iw)-wsub(ixo^s,iw))
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer, dimension(3, 3) kr
Kronecker delta tensor.
subroutine, public mhd_hllc_init()
subroutine mhd_get_wcd(wLC, wRC, whll, fRC, fLC, Fhll, patchf, lambdaCD, cmin, cmax, ixIL, ixOL, idim, f)
subroutine mhd_get_lcd(wLC, wRC, fLC, fRC, cmin, cmax, idim, ixIL, ixOL, whll, Fhll, lambdaCD, patchf)
subroutine mhd_diffuse_hllcd(ixIL, ixOL, idim, wLC, wRC, fLC, fRC, patchf)
Magneto-hydrodynamics module.
integer, dimension(:), allocatable, public, protected mom
Indices of the momentum density.
integer, public, protected rho_
Index of the density (in the w array)
procedure(sub_get_wcd), pointer phys_get_wcd
procedure(sub_get_lcd), pointer phys_get_lcd
procedure(sub_diffuse_hllcd), pointer phys_diffuse_hllcd