MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
mod_usr_template.t
Go to the documentation of this file.
1 !> This is a template for a new user problem of mhd
2 module mod_usr
3 
4  ! Include a physics module: mod_rho, mod_hd, mod_mhd ...
5  use mod_mhd
6 
7  implicit none
8 
9  ! Custom variables can be defined here
10  ! ...
11 
12 contains
13 
14  !> This routine should set user methods, and activate the physics module
15  subroutine usr_init()
16 
17  ! Choose coordinate system as 2D Cartesian with three components for vectors
18  call set_coordinate_system("Cartesian_2.5D")
19 
20  ! A routine for initial conditions is always required
21  usr_init_one_grid => initial_conditions
22 
23  ! Specify other user routines, for a list see mod_usr_methods.t
24  ! ...
25 
26  ! Choose independent normalization units if using dimensionless variables.
27  unit_length = 1.d9 ! cm
28  unit_temperature = 1.d6 ! K
29  unit_numberdensity = 1.d9 ! cm^-3
30 
31  ! Active the physics module: rho_activate(), hd_activate(), mhd_activate()
32  call mhd_activate()
33 
34  end subroutine usr_init
35 
36  !> A routine for specifying initial conditions
37  subroutine initial_conditions(ixI^L, ixO^L, w, x)
38 
39  integer, intent(in) :: ixI^L, ixO^L
40  double precision, intent(in) :: x(ixi^s,1:ndim)
41  double precision, intent(inout) :: w(ixi^s,1:nw)
42 
43  ! index for do loop and cell faces
44  integer :: idir, ixC^L
45 
46  ! Set initial values for w
47  ! 1.d0 and 0.d0 are just examples should be adjusted to user's problem
48 
49  ! cell-center density
50  w(ixo^s, rho_) = 1.d0
51  ! cell-center velocity
52  w(ixo^s, mom(1)) = 0.d0
53  w(ixo^s, mom(2)) = 0.d0
54  ! cell-center pressure
55  w(ixo^s, e_) = 1.d0
56  if(stagger_grid) then
57  ! set cell-face magnetic field B using CT method for divB control
58 
59  ! set B from vector potential (zero divB guaranteed) given
60  ! usr_init_vector_potential pointing to a subroutine to set vector potential
61  call b_from_vector_potential(ixgs^ll,ixi^l,ixo^l,block%ws,x)
62 
63  ! or directly set cell-face B (divB maybe non-zero) as following:
64  do idir=1,ndim
65  ixcmin^d=iximin^d;
66  ixcmax^d=iximax^d-kr(idir,^d);
67  ! cell-face B_idir
68  block%ws(ixc^s,idir)=1.d0
69  end do
70  ! update cell-center B from cell-face B
71  call mhd_face_to_center(ixo^l,block)
72  else
73  ! cell-center magnetic field
74  w(ixo^s,mag(1)) = 1.d0
75  w(ixo^s,mag(2)) = 1.d0
76  end if
77 
78  ! convert primitive variables to conservative variables
79  call mhd_to_conserved(ixi^l,ixo^l,w,x)
80 
81  end subroutine initial_conditions
82 
83  ! Extra routines can be placed here
84  ! ...
85 
86 end module mod_usr
Definition: mod_mhd.t:1
subroutine mhd_activate()
Definition: mod_mhd.t:15
This is a template for a new user problem of mhd.
subroutine initial_conditions(ixIL, ixOL, w, x)
A routine for specifying initial conditions.
subroutine usr_init()
This routine should set user methods, and activate the physics module.