MPI-AMRVAC  3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
mod_functions_bfield.t
Go to the documentation of this file.
2 
3  implicit none
4  private
5 
6 
7  public :: get_divb
8 
9  !> Indices of the magnetic field
10  integer, allocatable, public :: mag(:)
11 
12 contains
13 
14 
15  !> Calculate div B within ixO
16  subroutine get_divb(w,ixI^L,ixO^L,divb, fourthorder)
18  use mod_geometry
19 
20  integer, intent(in) :: ixi^l, ixo^l
21  double precision, intent(in) :: w(ixi^s,1:nw)
22  double precision, intent(inout) :: divb(ixi^s)
23  logical, intent(in), optional :: fourthorder
24 
25  integer :: ixc^l, idir
26 
27  if(stagger_grid) then
28  divb(ixo^s)=0.d0
29  do idir=1,ndim
30  ixc^l=ixo^l-kr(idir,^d);
31  divb(ixo^s)=divb(ixo^s)+block%ws(ixo^s,idir)*block%surfaceC(ixo^s,idir)-&
32  block%ws(ixc^s,idir)*block%surfaceC(ixc^s,idir)
33  end do
34  divb(ixo^s)=divb(ixo^s)/block%dvolume(ixo^s)
35  else
36  select case(typediv)
37  case("central")
38  call divvector(w(ixi^s,mag(1:ndir)),ixi^l,ixo^l,divb,fourthorder)
39  case("limited")
40  call divvectors(w(ixi^s,mag(1:ndir)),ixi^l,ixo^l,divb)
41  end select
42  end if
43 
44  end subroutine get_divb
45 
46 
47 end module mod_functions_bfield
integer, dimension(:), allocatable, public mag
Indices of the magnetic field.
subroutine, public get_divb(w, ixIL, ixOL, divb, fourthorder)
Calculate div B within ixO.
Module with geometry-related routines (e.g., divergence, curl)
Definition: mod_geometry.t:2
subroutine divvectors(qvec, ixIL, ixOL, divq)
Calculate divergence of a vector qvec within ixL using limited extrapolation to cell edges.
Definition: mod_geometry.t:608
subroutine divvector(qvec, ixIL, ixOL, divq, fourthorder, sixthorder)
Calculate divergence of a vector qvec within ixL.
Definition: mod_geometry.t:516
This module contains definitions of global parameters and variables and some generic functions/subrou...
type(state), pointer block
Block pointer for using one block and its previous state.
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.
character(len=std_len) typediv
integer, dimension(:), allocatable, parameter d
integer ndir
Number of spatial dimensions (components) for vector variables.