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