MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
amrini.t
Go to the documentation of this file.
1 !> Generate and initialize all grids at the coarsest level (level one)
2 subroutine initlevelone
5 
6  integer :: iigrid, igrid{#IFDEF EVOLVINGBOUNDARY , Morton_no}
7 
8  levmin=1
9  levmax=1
10 
11  call init_forest_root
12 
13  call getigrids
15 
16  ! fill solution space of all root grids
17  do iigrid=1,igridstail; igrid=igrids(iigrid);
18  call alloc_node(igrid)
19  ! in case gradient routine used in initial condition, ensure geometry known
20  call initial_condition(igrid)
21  end do
22  {#IFDEF EVOLVINGBOUNDARY
23  ! mark physical-boundary blocks on space-filling curve
24  do morton_no=morton_start(mype),morton_stop(mype)
25  igrid=sfc_to_igrid(morton_no)
26  if (phyboundblock(igrid)) sfc_phybound(morton_no)=1
27  end do
28  call mpi_allreduce(mpi_in_place,sfc_phybound,nleafs,mpi_integer,&
29  mpi_sum,icomm,ierrmpi)
30  }
31 
32  ! update ghost cells
33  call getbc(global_time,0.d0,ps,1,nwflux+nwaux)
34 
35 end subroutine initlevelone
36 
37 !> fill in initial condition
38 subroutine initial_condition(igrid)
39  ! Need only to set the mesh values (can leave ghost cells untouched)
42 
43  integer, intent(in) :: igrid
44 
45  ps(igrid)%w(ixg^t,1:nw)=zero
46 
47  saveigrid=igrid
48  ! in case gradient routine used in initial condition, ensure geometry known
49  block=>ps(igrid)
50  ^d&dxlevel(^d)=rnode(rpdx^d_,igrid);
53 
54  if (.not. associated(usr_init_one_grid)) then
55  call mpistop("usr_init_one_grid not defined")
56  else
57  call usr_init_one_grid(ixg^ll,ixm^ll,ps(igrid)%w,ps(igrid)%x)
58  end if
59 
60 end subroutine initial_condition
61 
62 !> modify initial condition
63 subroutine modify_ic
66 
67  integer :: iigrid, igrid
68 
69  do iigrid=1,igridstail; igrid=igrids(iigrid);
70  saveigrid=igrid
71  block=>ps(igrid)
72  ^d&dxlevel(^d)=rnode(rpdx^d_,igrid);
75 
76  if (.not. associated(usr_init_one_grid)) then
77  call mpistop("usr_init_one_grid not defined")
78  else
79  call usr_init_one_grid(ixg^ll,ixm^ll,ps(igrid)%w,ps(igrid)%x)
80  end if
81  end do
82 
83 end subroutine modify_ic
84 
85 {^nooned
86 !> improve initial condition after initialization
87 subroutine improve_initial_condition()
89  use mod_usr_methods
92  use mod_mhd_phys
94 
95  logical :: active
96 
97  if(associated(usr_improve_initial_condition)) then
99  else if(stagger_grid) then
100  if(associated(usr_init_vector_potential)) then
101  ! re-calculate magnetic field from the vector potential in a
102  ! completely divergency free way for AMR mesh in 3D
103  if(levmax>levmin.and.ndim==3) call recalculateb
104  end if
105  if(slab_uniform.and.clean_initial_divb) then
106  ! Project out the divB using multigrid poisson solver
107  ! if not initialised from vector potential
108  if(.not.use_multigrid) call mg_setup_multigrid()
109  call mhd_clean_divb_multigrid(global_time,0.d0,active)
110  call getbc(global_time,0.d0,ps,1,nwflux+nwaux)
111  end if
112  end if
113 
114 end subroutine improve_initial_condition
115 }
This module contains definitions of global parameters and variables and some generic functions/subrou...
subroutine, public mhd_clean_divb_multigrid(qdt, qt, active)
update ghost cells of all blocks including physical boundaries
logical slab_uniform
uniform Cartesian geometry or not (stretched Cartesian)
integer, dimension(:), allocatable type_gradient_limiter
Type of slope limiter used for computing gradients or divergences, when typegrad or typediv are set t...
integer, parameter plevel_
subroutine recalculateb
re-calculate the magnetic field from the vector potential in a completely divergency free way ...
integer, dimension(:), allocatable type_limiter
Type of slope limiter used for reconstructing variables on cell edges.
subroutine build_connectivity
Definition: connectivity.t:45
subroutine improve_initial_condition()
improve initial condition after initialization
Definition: amrini.t:88
Module to couple the octree-mg library to AMRVAC. This file uses the VACPP preprocessor, but its use is kept to a minimum.
subroutine init_forest_root
build root forest
Definition: forest.t:3
procedure(p_no_args), pointer usr_improve_initial_condition
Module with all the methods that users can customize in AMRVAC.
procedure(init_one_grid), pointer usr_init_one_grid
Initialize earch grid block data.
logical stagger_grid
True for using stagger grid.
integer ierrmpi
A global MPI error return code.
integer ixm
the mesh range (within a block with ghost cells)
subroutine getbc(time, qdt, psb, nwstart, nwbc, req_diag)
do update ghost cells of all blocks including physical boundaries
subroutine alloc_node(igrid)
allocate arrays on igrid node
integer, dimension(:), allocatable, parameter d
subroutine getigrids
Definition: connectivity.t:26
integer mype
The rank of the current MPI task.
double precision global_time
The global simulation time.
procedure(init_vector_potential), pointer usr_init_vector_potential
subroutine mpistop(message)
Exit MPI-AMRVAC with an error message.
Definition: comm_lib.t:198
subroutine modify_ic
modify initial condition
Definition: amrini.t:64
subroutine mg_setup_multigrid()
Setup multigrid for usage.
double precision, dimension(ndim) dxlevel
subroutine initial_condition(igrid)
fill in initial condition
Definition: amrini.t:39
integer, parameter ndim
Number of spatial dimensions for grid variables.
logical use_multigrid
Use multigrid (only available in 2D and 3D)
integer icomm
The MPI communicator.
integer, dimension(:,:), allocatable node
subroutine initlevelone
Generate and initialize all grids at the coarsest level (level one)
Definition: amrini.t:3
logical, dimension(:), allocatable phyboundblock
True if a block has any physical boundary.
Magneto-hydrodynamics module.
Definition: mod_mhd_phys.t:2
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
logical, public clean_initial_divb
clean initial divB
Definition: mod_mhd_phys.t:136