mod_afld Module Reference

Module for including anisotropic flux limited diffusion (AFLD)-approximation in Radiation-hydrodynamics simulations using mod_rhd Based on Turner and stone 2001. See [1]Moens, N., Sundqvist, J. O., El Mellah, I., Poniatowski, L., Teunissen, J., and Keppens, R., “Radiation-hydrodynamics with MPI-AMRVAC . Flux-limited diffusion”, Astronomy and Astrophysics, vol. 657, 2022. doi:10.1051/0004-6361/202141023. For more information. More...

Functions/Subroutines
subroutine	afld_params_read (files)
	public methods these are called in mod_rhd_phys
subroutine, public	afld_init (he_abundance, rhd_radiation_diffusion, afld_gamma)
	Initialising FLD-module: Read opacities Initialise Multigrid adimensionalise kappa Add extra variables to w-array, flux, kappa, eddington Tensor Lambda and R ...
subroutine, public	get_afld_rad_force (qdt, ixil, ixol, wct, w, x, energy, qsourcesplit, active)
	w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO This subroutine handles the radiation force
subroutine, public	afld_radforce_get_dt (w, ixil, ixol, dtnew, dxd, x)
subroutine, public	get_afld_energy_interact (qdt, ixil, ixol, wct, w, x, energy, qsourcesplit, active)
	w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO This subroutine handles the energy exchange between gas and radiation
subroutine, public	afld_get_opacity (w, x, ixil, ixol, fld_kappa)
	Sets the opacity in the w-array by calling mod_opal_opacity.
subroutine, public	afld_get_fluxlimiter (w, x, ixil, ixol, fld_lambda, fld_r, nth)
	Calculate fld flux limiter This subroutine calculates flux limiter lambda using the prescription stored in fld_fluxlimiter. It also calculates the ratio of radiation scaleheight and mean free path.
subroutine, public	afld_get_radflux (w, x, ixil, ixol, rad_flux, nth)
	Calculate Radiation Flux Calculate the Flux using the fld closure relation F = -c*lambda/(kappa*rho) *grad E.
subroutine	afld_get_eddington (w, x, ixil, ixol, eddington_tensor, nth)
	Calculate Eddington-tensor Stores Eddington-tensor in w-array.
subroutine, public	afld_get_radpress (w, x, ixil, ixol, rad_pressure, nth)
	Calculate Radiation Pressure Returns Radiation Pressure as tensor.
subroutine	diffuse_e_rad_mg (dtfactor, qdt, qtc, psa, psb)
	Calling all subroutines to perform the multigrid method Communicates rad_e and diff_coeff to multigrid library.
subroutine	evaluate_e_rad_mg (qtc, psa)
	inplace update of psa==>F_im(psa)
subroutine	put_diffterm_onegrid (ixil, ixol, w)
	inplace update of psa==>F_im(psa)
subroutine	afld_get_diffcoef_central (w, wct, x, ixil, ixol)
	Calculates cell-centered diffusion coefficient to be used in multigrid.
subroutine	update_diffcoeff (psa)
subroutine	afld_smooth_diffcoef (w, ixil, ixol, idir)
	Use running average on Diffusion coefficient.
subroutine	energy_interaction (w, wct, x, ixil, ixol)
	This subroutine calculates the radiation heating, radiation cooling and photon tiring using an implicit scheme. These sourceterms are applied using the root-finding of a 4th order polynomial This routine loops over every cell in the domain and computes the coefficients of the polynomials in every cell.
subroutine	bisection_method (e_gas, e_rad, c0, c1)
	Find the root of the 4th degree polynomial using the bisection method.
subroutine	newton_method (e_gas, e_rad, c0, c1)
	Find the root of the 4th degree polynomial using the Newton-Ralphson method.
subroutine	halley_method (e_gas, e_rad, c0, c1)
	Find the root of the 4th degree polynomial using the Halley method.
double precision function	polynomial_bisection (e_gas, c0, c1)
	Evaluate polynomial at argument e_gas.
double precision function	dpolynomial_bisection (e_gas, c0, c1)
	Evaluate first derivative of polynomial at argument e_gas.
double precision function	ddpolynomial_bisection (e_gas, c0, c1)
	Evaluate second derivative of polynomial at argument e_gas.

Variables
logical	fld_eint_split = .false.
	source split for energy interact and radforce:
logical	fld_radforce_split = .false.
double precision, public	fld_kappa0 = 0.34d0
	Opacity per unit of unit_density.
double precision, public	afld_mu = 0.6d0
	mean particle mass
double precision, public	fld_bisect_tol = 1.d-4
	Tolerance for bisection method for Energy sourceterms This is a percentage of the minimum of gas- and radiation energy.
double precision, public	fld_diff_tol = 1.d-4
	Tolerance for adi method for radiative Energy diffusion.
double precision, public	diff_crit
	Number for splitting the diffusion module.
character(len=8), dimension(:), allocatable	fld_opacity_law
	Use constant Opacity?
character(len=50)	fld_opal_table = 'Y09800'
character(len=16)	fld_fluxlimiter = 'Pomraning'
	Diffusion limit lambda = 0.33.
integer, dimension(:), allocatable	i_diff_mg
character(len=8)	afld_diff_scheme = 'mg'
	Which method to solve diffusion part.
character(len=8)	fld_interaction_method = 'Halley'
	Which method to find the root for the energy interaction polynomial.
logical	fld_diff_testcase = .false.
	Set Diffusion coefficient to unity.
logical	diff_coef_filter = .false.
	Take running average for Diffusion coefficient.
integer	size_d_filter = 1
logical	flux_lim_filter = .false.
	Take a running average over the fluxlimiter.
integer	size_l_filter = 1
logical	lineforce_opacities = .false.
	Use or dont use lineforce opacities.
logical	diffcrash_resume = .true.
	Resume run when multigrid returns error.
double precision	dt_diff = 0.d0
	running timestep for diffusion solver, initialised as zero

Detailed Description

Module for including anisotropic flux limited diffusion (AFLD)-approximation in Radiation-hydrodynamics simulations using mod_rhd Based on Turner and stone 2001. See [1]Moens, N., Sundqvist, J. O., El Mellah, I., Poniatowski, L., Teunissen, J., and Keppens, R., “Radiation-hydrodynamics with MPI-AMRVAC . Flux-limited diffusion”, Astronomy and Astrophysics, vol. 657, 2022. doi:10.1051/0004-6361/202141023. For more information.

Function/Subroutine Documentation

afld_get_diffcoef_central()

subroutine mod_afld::afld_get_diffcoef_central	(	double precision, dimension(ixi^s, 1:nw), intent(inout)	w,
		double precision, dimension(ixi^s, 1:nw), intent(in)	wct,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l
	)

Calculates cell-centered diffusion coefficient to be used in multigrid.

calculate diffusion coefficient

Definition at line 683 of file mod_afld.t.

Here is the call graph for this function:

afld_get_eddington()

subroutine mod_afld::afld_get_eddington	(	double precision, dimension(ixi^s, 1:nw), intent(in)	w,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixo^s,1:ndim,1:ndim), intent(out)	eddington_tensor,
		integer, intent(in)	nth
	)

Calculate Eddington-tensor Stores Eddington-tensor in w-array.

Calculate R everywhere |grad E|/(rho kappa E)

Calculate radiation pressure P = (lambda + lambda^2 R^2)*E

Definition at line 463 of file mod_afld.t.

Here is the call graph for this function:

afld_get_fluxlimiter()

subroutine, public mod_afld::afld_get_fluxlimiter	(	double precision, dimension(ixi^s, 1:nw), intent(in)	w,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixo^s,1:ndim), intent(out)	fld_lambda,
		double precision, dimension(ixo^s,1:ndim), intent(out)	fld_r,
		integer, intent(in)	nth
	)

Calculate fld flux limiter This subroutine calculates flux limiter lambda using the prescription stored in fld_fluxlimiter. It also calculates the ratio of radiation scaleheight and mean free path.

Calculate R everywhere |grad E|/(rho kappa E)

Calculate the flux limiter, lambda

Calculate R everywhere |grad E|/(rho kappa E)

Calculate the flux limiter, lambda Levermore and Pomraning: lambda = (2 + R)/(6 + 3R + R^2)

Calculate R everywhere |grad E|/(rho kappa E)

Calculate the flux limiter, lambda Minerbo:

Definition at line 337 of file mod_afld.t.

Here is the call graph for this function:

afld_get_opacity()

subroutine, public mod_afld::afld_get_opacity	(	double precision, dimension(ixi^s, 1:nw), intent(in)	w,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixo^s,1:ndim), intent(out)	fld_kappa
	)

Sets the opacity in the w-array by calling mod_opal_opacity.

Take lower value of rho in domain

Opacity bump

Take lower value of rho in domain

Hard limit on kappa

Definition at line 265 of file mod_afld.t.

Here is the call graph for this function:

afld_get_radflux()

subroutine, public mod_afld::afld_get_radflux	(	double precision, dimension(ixi^s, 1:nw), intent(in)	w,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixo^s, 1:ndim), intent(out)	rad_flux,
		integer, intent(in)	nth
	)

Calculate Radiation Flux Calculate the Flux using the fld closure relation F = -c*lambda/(kappa*rho) *grad E.

Definition at line 438 of file mod_afld.t.

Here is the call graph for this function:

afld_get_radpress()

subroutine, public mod_afld::afld_get_radpress	(	double precision, dimension(ixi^s, 1:nw), intent(in)	w,
		double precision, dimension(ixi^s, 1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixo^s,1:ndim,1:ndim), intent(out)	rad_pressure,
		integer, intent(in)	nth
	)

Calculate Radiation Pressure Returns Radiation Pressure as tensor.

Definition at line 517 of file mod_afld.t.

Here is the call graph for this function:

afld_init()

subroutine, public mod_afld::afld_init	(	double precision, intent(in)	he_abundance,
		logical, intent(in)	rhd_radiation_diffusion,
		double precision, intent(in)	afld_gamma
	)

Initialising FLD-module: Read opacities Initialise Multigrid adimensionalise kappa Add extra variables to w-array, flux, kappa, eddington Tensor Lambda and R ...

allocate boundary conditions

read par files

Need mean molecular weight

Definition at line 92 of file mod_afld.t.

Here is the call graph for this function:

afld_params_read()

subroutine mod_afld::afld_params_read

(

character(len=*), dimension(:), intent(in)

files

)

public methods these are called in mod_rhd_phys

Reading in fld-list parameters from .par file

Definition at line 66 of file mod_afld.t.

afld_radforce_get_dt()

subroutine, public mod_afld::afld_radforce_get_dt	(	double precision, dimension(ixi^s,1:nw), intent(in)	w,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, intent(inout)	dtnew,
		double precision, intent(in)	dx,
		double precision, intent(in)	d,
		double precision, dimension(ixi^s,1:ndim), intent(in)	x
	)

Definition at line 216 of file mod_afld.t.

Here is the call graph for this function:

afld_smooth_diffcoef()

subroutine mod_afld::afld_smooth_diffcoef	(	double precision, dimension(ixi^s, 1:nw), intent(inout)	w,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		integer, intent(in)	idir
	)

Use running average on Diffusion coefficient.

Definition at line 733 of file mod_afld.t.

bisection_method()

subroutine mod_afld::bisection_method	(	double precision, intent(inout)	e_gas,
		double precision, intent(in)	e_rad,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Find the root of the 4th degree polynomial using the bisection method.

Definition at line 843 of file mod_afld.t.

Here is the call graph for this function:

ddpolynomial_bisection()

double precision function mod_afld::ddpolynomial_bisection	(	double precision, intent(in)	e_gas,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Evaluate second derivative of polynomial at argument e_gas.

Definition at line 985 of file mod_afld.t.

diffuse_e_rad_mg()

subroutine mod_afld::diffuse_e_rad_mg	(	double precision, intent(in)	dtfactor,
		double precision, intent(in)	qdt,
		double precision, intent(in)	qtc,
		type(state), dimension(max_blocks), target	psa,
		type(state), dimension(max_blocks), target	psb
	)

Calling all subroutines to perform the multigrid method Communicates rad_e and diff_coeff to multigrid library.

Set diffusion timestep, add previous timestep if mg did not converge:

replace call set_rhs(mg, -1/dt, 0.0_dp)

Reset dt_diff when diffusion worked out

Definition at line 536 of file mod_afld.t.

Here is the call graph for this function:

dpolynomial_bisection()

double precision function mod_afld::dpolynomial_bisection	(	double precision, intent(in)	e_gas,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Evaluate first derivative of polynomial at argument e_gas.

Definition at line 975 of file mod_afld.t.

energy_interaction()

subroutine mod_afld::energy_interaction	(	double precision, dimension(ixi^s,1:nw), intent(inout)	w,
		double precision, dimension(ixi^s,1:nw), intent(in)	wct,
		double precision, dimension(ixi^s,1:ndim), intent(in)	x,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l
	)

This subroutine calculates the radiation heating, radiation cooling and photon tiring using an implicit scheme. These sourceterms are applied using the root-finding of a 4th order polynomial This routine loops over every cell in the domain and computes the coefficients of the polynomials in every cell.

e_gas is the INTERNAL ENERGY without KINETIC ENERGY

Calculate coefficients for polynomial

Loop over every cell for rootfinding method

advance E_rad

Update gas-energy in w, internal + kinetic

Beginning of module substracted wCT Ekin So now add wCT Ekin

Update rad-energy in w

Definition at line 764 of file mod_afld.t.

Here is the call graph for this function:

evaluate_e_rad_mg()

subroutine mod_afld::evaluate_e_rad_mg	(	double precision, intent(in)	qtc,
		type(state), dimension(max_blocks), target	psa
	)

inplace update of psa==>F_im(psa)

Definition at line 640 of file mod_afld.t.

Here is the call graph for this function:

get_afld_energy_interact()

subroutine, public mod_afld::get_afld_energy_interact	(	double precision, intent(in)	qdt,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixi^s,1:nw), intent(in)	wct,
		double precision, dimension(ixi^s,1:nw), intent(inout)	w,
		double precision, dimension(ixi^s,1:ndim), intent(in)	x,
		logical, intent(in)	energy,
		logical, intent(in)	qsourcesplit,
		logical, intent(inout)	active
	)

w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO This subroutine handles the energy exchange between gas and radiation

Calculate and add sourceterms

Add energy sourceterms

Definition at line 242 of file mod_afld.t.

Here is the call graph for this function:

get_afld_rad_force()

subroutine, public mod_afld::get_afld_rad_force	(	double precision, intent(in)	qdt,
		integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixi^s,1:nw), intent(in)	wct,
		double precision, dimension(ixi^s,1:nw), intent(inout)	w,
		double precision, dimension(ixi^s,1:ndim), intent(in)	x,
		logical, intent(in)	energy,
		logical, intent(in)	qsourcesplit,
		logical, intent(inout)	active
	)

w[iw]=w[iw]+qdt*S[wCT,qtC,x] where S is the source based on wCT within ixO This subroutine handles the radiation force

Calculate and add sourceterms

Radiation force = kappa*rho/c *Flux = lambda gradE

Momentum equation source term

Energy equation source term (kinetic energy)

Photon tiring calculate tensor div_v !$OMP PARALLEL DO

!$OMP END PARALLEL DO

VARIABLE NAMES DIV ARE ACTUALLY GRADIENTS

eq 34 Turner and stone (Only 2D)

Definition at line 139 of file mod_afld.t.

Here is the call graph for this function:

halley_method()

subroutine mod_afld::halley_method	(	double precision, intent(inout)	e_gas,
		double precision, intent(in)	e_rad,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Find the root of the 4th degree polynomial using the Halley method.

Compare error with dx = dx/dy dy

Definition at line 934 of file mod_afld.t.

Here is the call graph for this function:

newton_method()

subroutine mod_afld::newton_method	(	double precision, intent(inout)	e_gas,
		double precision, intent(in)	e_rad,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Find the root of the 4th degree polynomial using the Newton-Ralphson method.

Compare error with dx = dx/dy dy

Definition at line 905 of file mod_afld.t.

Here is the call graph for this function:

polynomial_bisection()

double precision function mod_afld::polynomial_bisection	(	double precision, intent(in)	e_gas,
		double precision, intent(in)	c0,
		double precision, intent(in)	c1
	)

Evaluate polynomial at argument e_gas.

Definition at line 965 of file mod_afld.t.

put_diffterm_onegrid()

subroutine mod_afld::put_diffterm_onegrid	(	integer, intent(in)	ixi,
		integer, intent(in)	l,
		integer, intent(in)	ixo,
			l,
		double precision, dimension(ixi^s, 1:nw), intent(inout)	w
	)

inplace update of psa==>F_im(psa)

Definition at line 662 of file mod_afld.t.

update_diffcoeff()

subroutine mod_afld::update_diffcoeff

(

type(state), dimension(max_blocks), target

psa

)

Definition at line 717 of file mod_afld.t.

Here is the call graph for this function:

Variable Documentation

afld_diff_scheme

character(len=8) mod_afld::afld_diff_scheme = 'mg'

Which method to solve diffusion part.

Definition at line 34 of file mod_afld.t.

afld_mu

double precision, public mod_afld::afld_mu = 0.6d0

mean particle mass

Definition at line 18 of file mod_afld.t.

diff_coef_filter

logical mod_afld::diff_coef_filter = .false.

Take running average for Diffusion coefficient.

Definition at line 40 of file mod_afld.t.

diff_crit

double precision, public mod_afld::diff_crit

Number for splitting the diffusion module.

Definition at line 25 of file mod_afld.t.

diffcrash_resume

logical mod_afld::diffcrash_resume = .true.

Resume run when multigrid returns error.

Definition at line 48 of file mod_afld.t.

dt_diff

double precision mod_afld::dt_diff = 0.d0

running timestep for diffusion solver, initialised as zero

Definition at line 52 of file mod_afld.t.

fld_bisect_tol

double precision, public mod_afld::fld_bisect_tol = 1.d-4

Tolerance for bisection method for Energy sourceterms This is a percentage of the minimum of gas- and radiation energy.

Definition at line 21 of file mod_afld.t.

fld_diff_testcase

logical mod_afld::fld_diff_testcase = .false.

Set Diffusion coefficient to unity.

Definition at line 38 of file mod_afld.t.

fld_diff_tol

double precision, public mod_afld::fld_diff_tol = 1.d-4

Tolerance for adi method for radiative Energy diffusion.

Definition at line 23 of file mod_afld.t.

fld_eint_split

logical mod_afld::fld_eint_split = .false.

source split for energy interact and radforce:

Definition at line 13 of file mod_afld.t.

fld_fluxlimiter

character(len=16) mod_afld::fld_fluxlimiter = 'Pomraning'

Diffusion limit lambda = 0.33.

Definition at line 30 of file mod_afld.t.

fld_interaction_method

character(len=8) mod_afld::fld_interaction_method = 'Halley'

Which method to find the root for the energy interaction polynomial.

Definition at line 36 of file mod_afld.t.

fld_kappa0

double precision, public mod_afld::fld_kappa0 = 0.34d0

Opacity per unit of unit_density.

Definition at line 16 of file mod_afld.t.

fld_opacity_law

character(len=8), dimension(:), allocatable mod_afld::fld_opacity_law

Use constant Opacity?

Definition at line 27 of file mod_afld.t.

fld_opal_table

character(len=50) mod_afld::fld_opal_table = 'Y09800'

Definition at line 28 of file mod_afld.t.

fld_radforce_split

logical mod_afld::fld_radforce_split = .false.

Definition at line 14 of file mod_afld.t.

flux_lim_filter

logical mod_afld::flux_lim_filter = .false.

Take a running average over the fluxlimiter.

Definition at line 43 of file mod_afld.t.

i_diff_mg

integer, dimension(:), allocatable mod_afld::i_diff_mg

Definition at line 32 of file mod_afld.t.

lineforce_opacities

logical mod_afld::lineforce_opacities = .false.

Use or dont use lineforce opacities.

Definition at line 46 of file mod_afld.t.

size_d_filter

integer mod_afld::size_d_filter = 1

Definition at line 41 of file mod_afld.t.

size_l_filter

integer mod_afld::size_l_filter = 1

Definition at line 44 of file mod_afld.t.