MPI-AMRVAC  3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
mod_dipole.t
Go to the documentation of this file.
1 !> sets up a magnetic dipole in a 3D cartesian box
2 ! parameters for dipole moment:
3 ! strength mu_dipole (factor 4pi included)
4 ! polar angle in degrees theta_dipole [0 to 90]
5 ! azimuthal angle in degrees phi_dipole [0 to 360]
6 ! center of dipole at x1_dip, x2_dip, x3_dip
7 module mod_dipole
8  implicit none
10  double precision :: mx, my, mz
11 
12 contains
13 
14 {^ifthreed
15  subroutine dipole(ixI^L,ixO^L,x,B)
17  use mod_comm_lib, only: mpistop
18 
19  integer, intent(in) :: ixI^L, ixO^L
20  double precision, intent(in) :: x(ixI^S,1:ndim)
21  double precision, intent(inout) :: B(ixI^S,1:ndir)
22 
23  double precision :: xsincos(ixI^S),xsinsin(ixI^S),xcos(ixI^S)
24  double precision :: rr0(ixI^S),mdotr(ixI^S)
25 
26  b=0.0d0
27  mx=mu_dipole*dsin(theta_dipole*dpi/180.0d0)*dcos(phi_dipole*dpi/180.d0)
28  my=mu_dipole*dsin(theta_dipole*dpi/180.0d0)*dsin(phi_dipole*dpi/180.0d0)
29  mz=mu_dipole*dcos(theta_dipole*dpi/180.0d0)
30  rr0(ixo^s)=dsqrt((x(ixo^s,1)-x1_dip)**2 &
31  +(x(ixo^s,2)-x2_dip)**2 &
32  +(x(ixo^s,3)-x3_dip)**2)
33  if(any(rr0(ixo^s)<smalldouble)) call mpistop("dipole center at cell center")
34  !where(rr0(ixO^S)<smalldouble)
35  ! rr0(ixO^S)=bigdouble
36  !endwhere
37  xsincos(ixo^s)=(x(ixo^s,1)-x1_dip)/rr0(ixo^s)
38  xsinsin(ixo^s)=(x(ixo^s,2)-x2_dip)/rr0(ixo^s)
39  xcos(ixo^s) =(x(ixo^s,3)-x3_dip)/rr0(ixo^s)
40  mdotr(ixo^s)=mx*xsincos(ixo^s)+my*xsinsin(ixo^s)+mz*xcos(ixo^s)
41  b(ixo^s,1)=(3.0d0*xsincos(ixo^s)*mdotr(ixo^s)-mx)/rr0(ixo^s)**3
42  b(ixo^s,2)=(3.0d0*xsinsin(ixo^s)*mdotr(ixo^s)-my)/rr0(ixo^s)**3
43  b(ixo^s,3)=(3.0d0*xcos(ixo^s) *mdotr(ixo^s)-mz)/rr0(ixo^s)**3
44 
45  end subroutine dipole
46 }
47 end module mod_dipole
subroutine, public mpistop(message)
Exit MPI-AMRVAC with an error message.
Definition: mod_comm_lib.t:208
sets up a magnetic dipole in a 3D cartesian box
Definition: mod_dipole.t:7
double precision x1_dip
Definition: mod_dipole.t:9
double precision mz
Definition: mod_dipole.t:10
double precision x3_dip
Definition: mod_dipole.t:9
double precision my
Definition: mod_dipole.t:10
double precision mx
Definition: mod_dipole.t:10
double precision phi_dipole
Definition: mod_dipole.t:9
double precision mu_dipole
Definition: mod_dipole.t:9
double precision x2_dip
Definition: mod_dipole.t:9
subroutine dipole(ixIL, ixOL, x, B)
Definition: mod_dipole.t:16
double precision theta_dipole
Definition: mod_dipole.t:9
This module contains definitions of global parameters and variables and some generic functions/subrou...