MPI-AMRVAC  2.0
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 
16  integer, intent(in):: ixI^L,ix^L
17  double precision, intent(inout) :: wB0(ixi^s,1:ndir)
18  double precision, intent(in) :: x(ixi^s,1:ndim)
19 
20  wb0(ix^s,1:ndir)=zero
21 
22  ! approximate cell-averaged B0 as cell-centered B0
23  select case (typeaxial)
24  case ("spherical")
25  {^nooned
26  if (dabs(bdip)>smalldouble) then
27  wb0(ix^s,1)=2.0d0*bdip*dcos(x(ix^s,2))/x(ix^s,1)**3
28  wb0(ix^s,2)=bdip*dsin(x(ix^s,2))/x(ix^s,1)**3
29  end if
30 
31  if (abs(bquad)>smalldouble) then
32  wb0(ix^s,1)=wb0(ix^s,1) &
33  +bquad*0.5d0*(1.0d0+3.0d0*dcos(2.0d0*x(ix^s,2)))/x(ix^s,1)**4
34  wb0(ix^s,2)=wb0(ix^s,2)+bquad*dsin(2.0d0*x(ix^s,2))/x(ix^s,1)**4
35  end if
36  if (abs(boct)>smalldouble) then
37  wb0(ix^s,1)=wb0(ix^s,1) &
38  +boct*(10.0d0*dcos(2.0d0*x(ix^s,2))-2.0d0) &
39  *dcos(x(ix^s,2))/x(ix^s,1)**5
40  wb0(ix^s,2)=wb0(ix^s,2) &
41  +boct*1.5d0*(3.0d0+5.0d0*dcos(2.0d0*x(ix^s,2))) &
42  *dsin(x(ix^s,2))/x(ix^s,1)**5
43  end if
44  }
45  end select
46 
47  if (associated(usr_set_b0)) call usr_set_b0(ixi^l,ix^l,x,wb0)
48 
49 end subroutine set_b0_cell
50 
51 subroutine set_j0_cell(igrid,wJ0,ixI^L,ix^L)
54  use mod_geometry
55 
56  integer, intent(in):: igrid,ixI^L,ix^L
57  double precision, intent(inout) :: wJ0(ixi^s,7-2*ndir:3)
58  integer :: idirmin0, idirmin
59 
60  if(associated(usr_set_j0)) then
61  call usr_set_j0(ixi^l,ix^l,ps(igrid)%x,wj0)
62  else
63  idirmin0 = 7-2*ndir
64  call curlvector(ps(igrid)%B0(:^d&,:,0),ixi^l,ix^l,wj0,idirmin,idirmin0,ndir)
65  end if
66 
67 end subroutine set_j0_cell
68 
69 subroutine set_b0_face(igrid,x,ixI^L,ix^L)
71 
72  integer, intent(in) :: igrid, ixI^L, ix^L
73  double precision, intent(in) :: x(ixi^s,1:ndim)
74 
75  double precision :: delx(ixi^s,1:ndim)
76  double precision :: xC(ixi^s,1:ndim),xshift^D
77  integer :: idims, ixC^L, ix, idims2
78 
79  if(slab)then
80  ^d&delx(ixi^s,^d)=rnode(rpdx^d_,igrid)\
81  else
82  ! for all non-cartesian and stretched coordinate(s)
83  delx(ixi^s,1:ndim)=ps(igrid)%dx(ixi^s,1:ndim)
84  endif
85 
86  do idims=1,ndim
87  ixcmin^d=ixmin^d-kr(^d,idims); ixcmax^d=ixmax^d;
88  ! always xshift=0 or 1/2
89  xshift^d=half*(one-kr(^d,idims));
90  do idims2=1,ndim
91  select case(idims2)
92  {case(^d)
93  do ix = ixc^lim^d
94  ! xshift=half: this is the cell center coordinate
95  ! xshift=0: this is the cell edge i+1/2 coordinate
96  xc(ix^d%ixC^s,^d)=x(ix^d%ixC^s,^d)+(half-xshift^d)*delx(ix^d%ixC^s,^d)
97  end do\}
98  end select
99  end do
100  call set_b0_cell(ps(igrid)%B0(:^d&,:,idims),xc,ixi^l,ixc^l)
101  end do
102 
103 end subroutine set_b0_face
104 
105 subroutine alloc_b0_grid(igrid)
107 
108  integer, intent(in) :: igrid
109 
110  if(.not. allocated(ps(igrid)%B0)) then
111  allocate(ps(igrid)%B0(ixg^t,1:ndir,0:ndim))
112  allocate(ps(igrid)%J0(ixg^t,7-2*ndir:3))
113  end if
114 
115 end subroutine alloc_b0_grid
116 
117 subroutine dealloc_b0_grid(igrid)
119 
120  integer, intent(in) :: igrid
121 
122  deallocate(ps(igrid)%B0)
123  deallocate(ps(igrid)%J0)
124 
125 end subroutine dealloc_b0_grid
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:70
subroutine set_j0_cell(igrid, wJ0, ixIL, ixL)
Definition: set_B0.t:52
Module with geometry-related routines (e.g., divergence, curl)
Definition: mod_geometry.t:2
integer ndir
Number of spatial dimensions (components) for vector variables.
character(len=std_len) typeaxial
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:485
integer, dimension(3, 3) kr
Kronecker delta tensor.
integer, dimension(:), allocatable, parameter d
subroutine dealloc_b0_grid(igrid)
Definition: set_B0.t:118
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.
logical slab
Cartesian geometry or not.
subroutine alloc_b0_grid(igrid)
Definition: set_B0.t:106
subroutine set_b0_cell(wB0, x, ixIL, ixL)
Definition: set_B0.t:13
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
procedure(set_j0), pointer usr_set_j0