MPI-AMRVAC 3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
Loading...
Searching...
No Matches
mod_b0.t
Go to the documentation of this file.
1module mod_b0
2 implicit none
3 private
4 public :: set_b0_grid
5
6contains
7
8 subroutine set_b0_grid(igrid)
10 integer, intent(in) :: igrid
11 integer :: ixcog^l
12
13 ixcogmin^d=1;
14 ixcogmax^d=(ixghi^d-2*nghostcells)/2+2*nghostcells;
15
16 call set_b0_cell(ps(igrid)%B0(ixg^t,:,0),ps(igrid)%x,ixg^ll,ixg^ll)
17 if(b0fieldalloccoarse) call set_b0_cell(psc(igrid)%B0(ixcog^s,:,0),psc(igrid)%x,ixcog^l,ixcog^l)
18 call set_j0_cell(igrid,ps(igrid)%J0,ixg^ll,ixm^ll^ladd1)
19 call set_b0_face(igrid,ps(igrid)%x,ixg^ll,ixm^ll)
20 end subroutine set_b0_grid
21
22 subroutine set_b0_cell(wB0,x,ixI^L,ix^L)
25 use mod_geometry
26
27 integer, intent(in):: ixi^l,ix^l
28 double precision, intent(inout) :: wb0(ixi^s,1:ndir)
29 double precision, intent(in) :: x(ixi^s,1:ndim)
30
31 wb0(ix^s,1:ndir)=zero
32
33 ! approximate cell-averaged B0 as cell-centered B0
34 select case (coordinate)
35 case (spherical)
36 {^nooned
37 if (dabs(bdip)>smalldouble) then
38 wb0(ix^s,1)=2.0d0*bdip*dcos(x(ix^s,2))/x(ix^s,1)**3
39 wb0(ix^s,2)=bdip*dsin(x(ix^s,2))/x(ix^s,1)**3
40 end if
41
42 if (abs(bquad)>smalldouble) then
43 wb0(ix^s,1)=wb0(ix^s,1) &
44 +bquad*0.5d0*(1.0d0+3.0d0*dcos(2.0d0*x(ix^s,2)))/x(ix^s,1)**4
45 wb0(ix^s,2)=wb0(ix^s,2)+bquad*dsin(2.0d0*x(ix^s,2))/x(ix^s,1)**4
46 end if
47 if (abs(boct)>smalldouble) then
48 wb0(ix^s,1)=wb0(ix^s,1) &
49 +boct*(10.0d0*dcos(2.0d0*x(ix^s,2))-2.0d0) &
50 *dcos(x(ix^s,2))/x(ix^s,1)**5
51 wb0(ix^s,2)=wb0(ix^s,2) &
52 +boct*1.5d0*(3.0d0+5.0d0*dcos(2.0d0*x(ix^s,2))) &
53 *dsin(x(ix^s,2))/x(ix^s,1)**5
54 end if
55 }
56 end select
57 if (associated(usr_set_b0)) call usr_set_b0(ixi^l,ix^l,x,wb0)
58 end subroutine set_b0_cell
59
60 subroutine set_j0_cell(igrid,wJ0,ixI^L,ix^L)
63 use mod_geometry
64
65 integer, intent(in):: igrid,ixi^l,ix^l
66 double precision, intent(inout) :: wj0(ixi^s,7-2*ndir:3)
67 integer :: idirmin0, idirmin
68
69 if(associated(usr_set_j0)) then
70 call usr_set_j0(ixi^l,ix^l,ps(igrid)%x,wj0)
71 else
72 idirmin0 = 7-2*ndir
73 call curlvector(ps(igrid)%B0(ixi^s,:,0),ixi^l,ix^l,wj0,idirmin,idirmin0,ndir)
74 end if
75 end subroutine set_j0_cell
76
77 subroutine set_b0_face(igrid,x,ixI^L,ixO^L)
79
80 integer, intent(in) :: igrid, ixi^l, ixo^l
81 double precision, intent(in) :: x(ixi^s,1:ndim)
82
83 double precision :: delx(ixi^s,1:ndim)
84 double precision :: xc(ixi^s,1:ndim),xshift^d
85 integer :: idims, ixc^l, hxo^l, ix, idims2
86
87 if(slab_uniform)then
88 ^d&delx(ixi^s,^d)=rnode(rpdx^d_,igrid)\
89 else
90 ! for all non-cartesian and stretched cartesian coordinates
91 delx(ixi^s,1:ndim)=ps(igrid)%dx(ixi^s,1:ndim)
92 endif
93
94 do idims=1,ndim
95 hxo^l=ixo^l-kr(idims,^d);
96 if(stagger_grid) then
97 ! ct needs all transverse cells
98 ixcmax^d=ixomax^d+nghostcells-nghostcells*kr(idims,^d); ixcmin^d=hxomin^d-nghostcells+nghostcells*kr(idims,^d);
99 else
100 ! ixC is centered index in the idims direction from ixOmin-1/2 to ixOmax+1/2
101 ixcmax^d=ixomax^d; ixcmin^d=hxomin^d;
102 end if
103 ! always xshift=0 or 1/2
104 xshift^d=half*(one-kr(^d,idims));
105 do idims2=1,ndim
106 select case(idims2)
107 {case(^d)
108 do ix = ixc^lim^d
109 ! xshift=half: this is the cell center coordinate
110 ! xshift=0: this is the cell edge i+1/2 coordinate
111 xc(ix^d%ixC^s,^d)=x(ix^d%ixC^s,^d)+(half-xshift^d)*delx(ix^d%ixC^s,^d)
112 end do\}
113 end select
114 end do
115 call set_b0_cell(ps(igrid)%B0(ixi^s,:,idims),xc,ixi^l,ixc^l)
116 end do
117 end subroutine set_b0_face
118end module mod_b0
119
Definition mod_b0.t:1
subroutine, public set_b0_grid(igrid)
Definition mod_b0.t:9
Module with geometry-related routines (e.g., divergence, curl)
Definition mod_geometry.t:2
integer coordinate
Definition mod_geometry.t:7
integer, parameter spherical
subroutine curlvector(qvec, ixil, ixol, curlvec, idirmin, idirmin0, ndir0, fourthorder)
Calculate curl of a vector qvec within ixL Options to employ standard second order CD evaluations use...
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer ixghi
Upper index of grid block arrays.
integer, dimension(3, 3) kr
Kronecker delta tensor.
integer, parameter ndim
Number of spatial dimensions for grid variables.
logical stagger_grid
True for using stagger grid.
double precision bdip
amplitude of background dipolar, quadrupolar, octupolar, user's field
double precision, dimension(:), allocatable, parameter d
integer ndir
Number of spatial dimensions (components) for vector variables.
integer ixm
the mesh range of a physical block without ghost cells
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
integer nghostcells
Number of ghost cells surrounding a grid.
logical slab_uniform
uniform Cartesian geometry or not (stretched Cartesian)
Module with all the methods that users can customize in AMRVAC.
procedure(set_j0), pointer usr_set_j0
procedure(set_b0), pointer usr_set_b0