! $Id$
!
! This module provide a way for users to specify custom
! (i.e. not in the standard Pencil Code) physics, diagnostics etc.
!
! The module provides a set of standard hooks into the Pencil-Code and
! currently allows the following customizations:
!
! Description | Relevant function call
! ---------------------------------------------------------------------------
! Special variable registration | register_special
! (pre parameter read) |
! Special variable initialization | initialize_special
! (post parameter read) |
! Special variable finalization | finalize_special
! (deallocation, etc.) |
! |
! Special initial condition | init_special
! this is called last so may be used to modify |
! the mvar variables declared by this module |
! or optionally modify any of the other f array |
! variables. The latter, however, should be |
! avoided where ever possible. |
! |
! Special term in the mass (density) equation | special_calc_density
! Special term in the momentum (hydro) equation | special_calc_hydro
! Special term in the energy equation | special_calc_energy
! Special term in the induction (magnetic) | special_calc_magnetic
! equation |
! |
! Special equation | dspecial_dt
! NOT IMPLEMENTED FULLY YET - HOOKS NOT PLACED INTO THE PENCIL-CODE
!
!** AUTOMATIC CPARAM.INC GENERATION ****************************
! Declare (for generation of special_dummies.inc) the number of f array
! variables and auxiliary variables added by this module
!
! CPARAM logical, parameter :: lspecial = .true.
!
! MVAR CONTRIBUTION 0
! MAUX CONTRIBUTION 0
!
!***************************************************************
!
! HOW TO USE THIS FILE
! --------------------
!
! Change the line above to
! lspecial = .true.
! to enable use of special hooks.
!
! The rest of this file may be used as a template for your own
! special module. Lines which are double commented are intended
! as examples of code. Simply fill out the prototypes for the
! features you want to use.
!
! Save the file with a meaningful name, eg. geo_kws.f90 and place
! it in the $PENCIL_HOME/src/special directory. This path has
! been created to allow users ot optionally check their contributions
! in to the Pencil-Code SVN repository. This may be useful if you
! are working on/using the additional physics with somebodyelse or
! may require some assistance from one of the main Pencil-Code team.
!
! To use your additional physics code edit the Makefile.local in
! the src directory under the run directory in which you wish to
! use your additional physics. Add a line with all the module
! selections to say something like:
!
! SPECIAL=special/geo_kws
!
! Where geo_kws it replaced by the filename of your new module
! upto and not including the .f90
!
module Special
!
use Cparam
use Cdata
use General, only: keep_compiler_quiet
use Messages, only: svn_id, fatal_error
!
implicit none
!
include '../special.h'
!
real :: gravz_boussinesq=0., gravz_boussinesq_lnrho=0., betaz_boussinesq=0.
!
namelist /special_run_pars/ gravz_boussinesq, gravz_boussinesq_lnrho, &
betaz_boussinesq
!
! Declare index of new variables in f array (if any).
!
!! integer :: ispecial=0
!! integer :: ispecaux=0
!
!! Diagnostic variables (needs to be consistent with reset list below).
!
!! integer :: idiag_POSSIBLEDIAGNOSTIC=0
!
contains
!***********************************************************************
subroutine register_special()
!
! Set up indices for variables in special modules.
!
! 6-oct-03/tony: coded
!
if (lroot) call svn_id( &
"$Id$")
!
!! call farray_register_pde('special',ispecial)
!! call farray_register_auxiliary('specaux',ispecaux)
!! call farray_register_auxiliary('specaux',ispecaux,communicated=.true.)
!
endsubroutine register_special
!***********************************************************************
subroutine initialize_special(f)
!
! Called after reading parameters, but before the time loop.
!
! 06-oct-03/tony: coded
!
real, dimension (mx,my,mz,mfarray) :: f
!
call keep_compiler_quiet(f)
!
endsubroutine initialize_special
!***********************************************************************
subroutine finalize_special(f)
!
! Called right before exiting.
!
! 14-aug-2011/Bourdin.KIS: coded
!
real, dimension (mx,my,mz,mfarray), intent(inout) :: f
!
call keep_compiler_quiet(f)
!
endsubroutine finalize_special
!***********************************************************************
subroutine init_special(f)
!
! initialise special condition; called from start.f90
! 06-oct-2003/tony: coded
!
real, dimension (mx,my,mz,mfarray) :: f
!
intent(inout) :: f
!!
!! SAMPLE IMPLEMENTATION
!!
!! select case (initspecial)
!! case ('nothing'); if (lroot) print*,'init_special: nothing'
!! case ('zero', '0'); f(:,:,:,iSPECIAL_VARIABLE_INDEX) = 0.
!! case default
!! call fatal_error("init_special: No such value for initspecial:" &
!! ,trim(initspecial))
!! endselect
!
call keep_compiler_quiet(f)
!
endsubroutine init_special
!***********************************************************************
subroutine pencil_criteria_special()
!
! All pencils that this special module depends on are specified here.
!
! 25-dec-14/axel: adapted from nospecial
!
if (lhydro.and.lentropy) then
lpenc_requested(i_cp)=.true.
lpenc_requested(i_uu)=.true.
if (gravz_boussinesq/=0.) lpenc_requested(i_ss)=.true.
if (gravz_boussinesq/=0.) lpenc_requested(i_cp1)=.true.
if (gravz_boussinesq_lnrho/=0.) lpenc_requested(i_lnrho)=.true.
endif
!
endsubroutine pencil_criteria_special
!***********************************************************************
subroutine pencil_interdep_special(lpencil_in)
!
! Interdependency among pencils provided by this module are specified here.
!
! 18-07-06/tony: coded
!
logical, dimension(npencils), intent(inout) :: lpencil_in
!
call keep_compiler_quiet(lpencil_in)
!
endsubroutine pencil_interdep_special
!***********************************************************************
subroutine calc_pencils_special(f,p)
!
! Calculate Special pencils.
! Most basic pencils should come first, as others may depend on them.
!
! 24-nov-04/tony: coded
!
real, dimension (mx,my,mz,mfarray) :: f
type (pencil_case) :: p
!
intent(in) :: f
intent(inout) :: p
!
call keep_compiler_quiet(f)
call keep_compiler_quiet(p)
!
endsubroutine calc_pencils_special
!***********************************************************************
subroutine dspecial_dt(f,df,p)
!
! calculate right hand side of ONE OR MORE extra coupled PDEs
! along the 'current' Pencil, i.e. f(l1:l2,m,n) where
! m,n are global variables looped over in equ.f90
!
! Due to the multi-step Runge Kutta timestepping used one MUST always
! add to the present contents of the df array. NEVER reset it to zero.
!
! Several precalculated Pencils of information are passed for
! efficiency.
!
! 06-oct-03/tony: coded
!
real, dimension (mx,my,mz,mfarray) :: f
real, dimension (mx,my,mz,mvar) :: df
type (pencil_case) :: p
!
intent(in) :: f,p
intent(inout) :: df
!
! Identify module and boundary conditions.
!
if (headtt.or.ldebug) print*,'dspecial_dt: SOLVE dspecial_dt'
!! if (headtt) call identify_bcs('special',ispecial)
!
!!
!! SAMPLE DIAGNOSTIC IMPLEMENTATION
!!
!! if (ldiagnos) then
!! if (idiag_SPECIAL_DIAGNOSTIC/=0) then
!! call sum_mn_name(MATHEMATICAL EXPRESSION,idiag_SPECIAL_DIAGNOSTIC)
!!! see also integrate_mn_name
!! endif
!! endif
!
call keep_compiler_quiet(f,df)
call keep_compiler_quiet(p)
!
endsubroutine dspecial_dt
!***********************************************************************
subroutine read_special_run_pars(iostat)
!
use File_io, only: parallel_unit
!
integer, intent(out) :: iostat
!
read(parallel_unit, NML=special_run_pars, IOSTAT=iostat)
!
endsubroutine read_special_run_pars
!***********************************************************************
subroutine write_special_run_pars(unit)
!
integer, intent(in) :: unit
!
write(unit, NML=special_run_pars)
!
endsubroutine write_special_run_pars
!***********************************************************************
subroutine rprint_special(lreset,lwrite)
!
! Reads and registers print parameters relevant to special.
!
! 06-oct-03/tony: coded
!
!! use FArrayManager, only: farray_index_append
!
!! integer :: iname
logical :: lreset,lwr
logical, optional :: lwrite
!
lwr = .false.
if (present(lwrite)) lwr=lwrite
!!!
!!! reset everything in case of reset
!!! (this needs to be consistent with what is defined above!)
!!!
if (lreset) then
!! idiag_SPECIAL_DIAGNOSTIC=0
endif
!!
!! do iname=1,nname
!! call parse_name(iname,cname(iname),cform(iname),&
!! 'NAMEOFSPECIALDIAGNOSTIC',idiag_SPECIAL_DIAGNOSTIC)
!! enddo
!!
!!! write column where which magnetic variable is stored
!! if (lwr) then
!! call farray_index_append('idiag_SPECIAL_DIAGNOSTIC',idiag_SPECIAL_DIAGNOSTIC)
!! endif
!!
endsubroutine rprint_special
!***********************************************************************
subroutine get_slices_special(f,slices)
!
! Write slices for animation of Special variables.
!
! 26-jun-06/tony: dummy
!
real, dimension (mx,my,mz,mfarray) :: f
type (slice_data) :: slices
!
call keep_compiler_quiet(f)
call keep_compiler_quiet(slices%ready)
!
endsubroutine get_slices_special
!***********************************************************************
subroutine special_calc_hydro(f,df,p)
!
! Added Boussinesq term, assuming z-direction in Cartesian case
!
! 25-dec-14/axel: adapted from nospecial
!
real, dimension (mx,my,mz,mfarray), intent(in) :: f
real, dimension (mx,my,mz,mvar), intent(inout) :: df
type (pencil_case), intent(in) :: p
!
! gravz_boussineq term
!
if (lhydro.and.lentropy) then
if (gravz_boussinesq_lnrho/=0.) then
df(l1:l2,m,n,iuz)=df(l1:l2,m,n,iuz)-p%lnrho*gravz_boussinesq_lnrho
elseif (gravz_boussinesq/=0.) then
df(l1:l2,m,n,iuz)=df(l1:l2,m,n,iuz)+p%ss*p%cp1*gravz_boussinesq
endif
endif
!
call keep_compiler_quiet(f)
!
endsubroutine special_calc_hydro
!***********************************************************************
subroutine special_calc_density(f,df,p)
!
! Calculate an additional 'special' term on the right hand side of the
! continuity equation.
!
! Some precalculated pencils of data are passed in for efficiency
! others may be calculated directly from the f array
!
! 06-oct-03/tony: coded
!
real, dimension (mx,my,mz,mfarray), intent(in) :: f
real, dimension (mx,my,mz,mvar), intent(inout) :: df
type (pencil_case), intent(in) :: p
!!
!! SAMPLE IMPLEMENTATION (remember one must ALWAYS add to df).
!!
!! df(l1:l2,m,n,ilnrho) = df(l1:l2,m,n,ilnrho) + SOME NEW TERM
!!
call keep_compiler_quiet(f,df)
call keep_compiler_quiet(p)
!
endsubroutine special_calc_density
!***********************************************************************
subroutine special_calc_dustdensity(f,df,p)
!
! Calculate an additional 'special' term on the right hand side of the
! continuity equation.
!
! Some precalculated pencils of data are passed in for efficiency
! others may be calculated directly from the f array
!
! 06-oct-03/tony: coded
!
real, dimension (mx,my,mz,mfarray), intent(in) :: f
real, dimension (mx,my,mz,mvar), intent(inout) :: df
type (pencil_case), intent(in) :: p
!!
!! SAMPLE IMPLEMENTATION (remember one must ALWAYS add to df).
!!
!! df(l1:l2,m,n,ilnrho) = df(l1:l2,m,n,ilnrho) + SOME NEW TERM
!!
call keep_compiler_quiet(f,df)
call keep_compiler_quiet(p)
!
endsubroutine special_calc_dustdensity
!***********************************************************************
subroutine special_calc_energy(f,df,p)
!
! Added Boussinesq term, assuming z-direction in Cartesian case
!
! 25-dec-14/axel: adapted from nospecial
!
real, dimension (mx,my,mz,mfarray), intent(in) :: f
real, dimension (mx,my,mz,mvar), intent(inout) :: df
type (pencil_case), intent(in) :: p
!
! betaz_boussinesq term
!
if (lhydro.and.lentropy) then
df(l1:l2,m,n,iss) = df(l1:l2,m,n,iss) + p%uu(:,3)*p%cp*betaz_boussinesq
endif
!
call keep_compiler_quiet(f)
!
endsubroutine special_calc_energy
!***********************************************************************
subroutine special_calc_magnetic(f,df,p)
!
! Calculate an additional 'special' term on the right hand side of the
! induction equation.
!
! Some precalculated pencils of data are passed in for efficiency
! others may be calculated directly from the f array.
!
! 06-oct-03/tony: coded
!
real, dimension (mx,my,mz,mfarray), intent(in) :: f
real, dimension (mx,my,mz,mvar), intent(inout) :: df
type (pencil_case), intent(in) :: p
!!
!! SAMPLE IMPLEMENTATION (remember one must ALWAYS add to df).
!!
!! df(l1:l2,m,n,iux) = df(l1:l2,m,n,iux) + SOME NEW TERM
!! df(l1:l2,m,n,iuy) = df(l1:l2,m,n,iuy) + SOME NEW TERM
!! df(l1:l2,m,n,iuz) = df(l1:l2,m,n,iuz) + SOME NEW TERM
!!
call keep_compiler_quiet(f,df)
call keep_compiler_quiet(p)
!
endsubroutine special_calc_magnetic
!***********************************************************************
subroutine special_calc_pscalar(f,df,p)
!
! Calculate an additional 'special' term on the right hand side of the
! passive scalar equation.
!
! Some precalculated pencils of data are passed in for efficiency
! others may be calculated directly from the f array.
!
! 15-jun-09/anders: coded
!
real, dimension (mx,my,mz,mfarray), intent(in) :: f
real, dimension (mx,my,mz,mvar), intent(inout) :: df
type (pencil_case), intent(in) :: p
!!
!! SAMPLE IMPLEMENTATION (remember one must ALWAYS add to df).
!!
!! df(l1:l2,m,n,ilncc) = df(l1:l2,m,n,ilncc) + SOME NEW TERM
!!
call keep_compiler_quiet(f,df)
call keep_compiler_quiet(p)
!
endsubroutine special_calc_pscalar
!***********************************************************************
subroutine special_calc_particles(fp)
!
! Called before the loop, in case some particle value is needed
! for the special density/hydro/magnetic/entropy.
!
! 20-nov-08/wlad: coded
!
real, dimension (:,:), intent(in) :: fp
!
call keep_compiler_quiet(fp)
!
endsubroutine special_calc_particles
!***********************************************************************
subroutine special_calc_chemistry(f,df,p)
!
! Calculate an additional 'special' term on the right hand side of the
! induction equation.
!
!
! 15-sep-10/natalia: coded
!
real, dimension (mx,my,mz,mfarray), intent(in) :: f
real, dimension (mx,my,mz,mvar), intent(inout) :: df
type (pencil_case), intent(in) :: p
!!
!! SAMPLE IMPLEMENTATION (remember one must ALWAYS add to df).
!!
!! df(l1:l2,m,n,iux) = df(l1:l2,m,n,iux) + SOME NEW TERM
!! df(l1:l2,m,n,iuy) = df(l1:l2,m,n,iuy) + SOME NEW TERM
!! df(l1:l2,m,n,iuz) = df(l1:l2,m,n,iuz) + SOME NEW TERM
!!
call keep_compiler_quiet(f,df)
call keep_compiler_quiet(p)
!
endsubroutine special_calc_chemistry
!***********************************************************************
subroutine special_before_boundary(f)
!
! Possibility to modify the f array before the boundaries are
! communicated.
!
! Some precalculated pencils of data are passed in for efficiency
! others may be calculated directly from the f array
!
! 06-jul-06/tony: coded
!
real, dimension (mx,my,mz,mfarray), intent(in) :: f
!
call keep_compiler_quiet(f)
!
endsubroutine special_before_boundary
!***********************************************************************
subroutine special_after_timestep(f,df,dt_,llast)
!
! Possibility to modify the f and df after df is updated.
! Used for the Fargo shift, for instance.
!
! 27-nov-08/wlad: coded
!
logical, intent(in) :: llast
real, dimension(mx,my,mz,mfarray), intent(inout) :: f
real, dimension(mx,my,mz,mvar), intent(inout) :: df
real, intent(in) :: dt_
!
call keep_compiler_quiet(f,df)
call keep_compiler_quiet(dt_)
call keep_compiler_quiet(llast)
!
endsubroutine special_after_timestep
!***********************************************************************
subroutine set_init_parameters(Ntot,BB0,dsize,init_distr,init_distr2)
!
! Possibility to modify the f and df after df is updated.
! Used for the Fargo shift, for instance.
!
! 27-nov-08/wlad: coded
!
real, dimension(mx,my,mz,mfarray) :: f
real, dimension(ndustspec) :: dsize,init_distr,init_distr2
real :: Ntot, BB0
!
call keep_compiler_quiet(f)
call keep_compiler_quiet(dsize,init_distr,init_distr2)
call keep_compiler_quiet(Ntot,BB0)
!
endsubroutine set_init_parameters
!***********************************************************************
endmodule Special