** Simple Disk Envelope (SiDE) ﬁt framework**

**4.6 Example of basic usage 87**

*# Set image parameters*

impar = [{'npix':512,'wav':[1100.,3000.],'sizeau':11000,'incl':67.}]

• **Model image options** In visdata and impar have to be deﬁne wav which is the
wavelength in microns, npix is the number of pixels on the rectangular images,
sizeau is the diameter of the image in au and incl is the inclination angle of the
source. In this framework the input original image size can be re-scaled to speed-up
the model image time computing. This was implemented in the previous analysis
of Chapter 3, where we required a large image size to model the observation, but
only a small percent from extended emission contributed to the modelled *uv-range.*

To avoid do the raytracing on pixels that doesn’t contribute much to the overall visibilitysizeau_extis used. The new image conserves the original image resolution (in terms of au/px), the old image size is speciﬁed insizeau_ext andsizeau covers the new image size. The original image is stored in the center of the new image, while the extended image size is ﬁlled out with zero values. It is recommended to check the data-sets withGALARIO, using get_image_size, before to implementing this change. In case of a typical model using for example the data-sets of Chapter 3, the computation time is dominated by the raytracing step. To fully resolve the models on the same scale as the observations, large synthetic images are needed. Typical runtime of the thermal radiation transport is 30 – 40 s, forGALARIO runtime is 6 – 10 s.

• **Model options** In the paramfile is possible to setup the stellar properties, dust
opacity and density proﬁles for diﬀerent components. The disk model is given by a
ﬂaring density proﬁle in which is possible to change parameters as the height of the
disk, ﬂaring index, density power law, inner disk radius among others. The envelope
proﬁle has four options: powerlaw, Tafalla2004, Ulrich1976andUlrich1976_2.

In Table 4.1 are the parameters that have to be setup for each envelope density pro-ﬁle. In the case of powerlaw and Tafalla2004density models the parameterr0Env refers to the radius at which rho0Envis deﬁned, while in the case of theUlrich1976 density proﬁle r0Env is the centrifugal radius.

The block of dust opacity contains the option to change the number of grain pop-ulations, ngpop. This version of the framework accepts up to three diﬀerent dust populations, one for the disk and two for the envelope. If Ulrich1976_2 is setup, Envcuthas to be deﬁned in theparamfile. Envcutis a limit in units of au where one dust population is in the inner part of Envcut and another dust population outside Envcut.

Additionally, if the model contains an outﬂow cavity then it is possible to setup dif-ferent cavity shapes with modeCav. The available cavity shapes are: Sheehan2017, edgecenter and gridcenter, see Fig. 4.2. Sheehan2017 is a outﬂow cavity model from Sheehan & Eisner (2017) characterized by:

*z* ą1 au`*r*^{ζ}*,* (4.2)

Table 4.1: Model setup parameters Parameter

Envelope model rho0Env r0Env* ^{a}* rTrunEnv prhoEnv ngpop

powerlaw 3 3 3 3 1

Tafalla2004 3 3 3 3 1

Ulrich1976 3 3 3 7 1

Ulrich1976_2* ^{b}* 3 3 3 7 2

(* ^{a}*)r0Env is the centrifugal radius in the Ulrich proﬁle. (

*) To setup 2 diﬀerent dust populations in the Ulrich proﬁle Envcut has to be deﬁne in units of au.*

^{b}where 1 au is the height above the midplane where the cavity starts and *ζ*
charac-terizes the shape and opening angle of the cavity. More details in Eisner (2012).

In the gridcenter option the cavity start at the center of the model grid (0,0,0), while in theedgecenter the cavity start at the outer radius of the disk following the equation:

arctan

ˆ*r*sin*θ*´*r*_{disk}*r*cos*θ*

˙

ă*θ*_{cav}*,* (4.3)

where *r*_{disk} is the outer radius of the disk in au, and *θ*_{cav} is the open angle onf the
cavity in degrees.

• **Fitting arguments** The next step is to deﬁne a dictionary that contains keyword
arguments for the logarithm of the posterior probability function lnpostfn(). An
example of kwargs is given by:

kwargs = {'dpc': 125., 'incl': 67., 'impar': impar, 'verbose': True, 'idisk':True, 'ienv':True, 'icav':False,

'cleanModel': True, 'binary': True, 'chi2_only':True, 'galario_check':False, 'time':True }

where dpc is the distace to the source in unit of parsec, impar collects the image parameters, verbose is an option to print a summary of the model parameters and the inclusion of the model proﬁles withidisk (disk proﬁle),ienv(envelope proﬁle), icav (outﬂow cavity) and islab (slab model). cleanModel is an option to delete the RADMC-3D model folder from disk after the posterior probability estimation, therefore the model ﬁles are not stored. If binary is True then RADMC3D will use binary format otherwise the ﬁle format is ascii. The use of binary ﬁles improves the computation speed and reduces disk space usage when models are kept.

The keyword chi2_onlycomputes and store the *χ*^{2} between the data and model. If
it is True the synthetic visibility is computed but not stored (a zero value array is
stored instead). The*χ*^{2} is still computed and stored. Set this to True when running
MCMC ﬁtting in order to improve speed.

Figure 4.2: Three examples of 2D density maps for diﬀerent shapes of outﬂow cavity. In gridcenter the outﬂow cavity starts at the center of the grid, edgecentercavity starts atrTrunEnv, the inner radius of the envelope density proﬁle and the Sheehan2017is a cavity described in Sheehan & Eisner (2017).

Finally, galario_check and time, return the Nyquist criterion for computing the
synthetic visibilities in the (u,*v) locations provided by GALARIO and the function*
run-time information, respectively.

The following step is to choose the ﬁtting parameters, the range of the space ex-ploration and the initial guess for the parameters. The ﬁtted parameter names are listed in parname. These can be parameters from the model setup ﬁle and some special parameters as: source distance, inclination, position angle and oﬀset in right ascension (RA) and declination (DEC). p_rangesare the limits for each parameter.

*# Choose fitting parameters and initial values*

parname = ['mdisk','rho0Env','gsmax_disk','gsmax_env']

p_ranges = [[-10., -2.], *# log disk dust mass [solar mass]*

[-23., -19.], *# log envelope dust density [g/cm**3]*

[-6., 0.], *# log disk grain size [cm]*

[-6., 0.]] *# log envelope grain size [cm]*

Additionally, the initial guess for the parameters is assigned to p0. The mode of exploration of the initial guess parameter can be deﬁned through the lists: p_form, p_formpriorand p_sigma.

*# initial guess for the parameters*
p0 = [-5, -20, -4., -4.]

p_form = ['log', 'lin', 'log', 'log']

p_formprior = ['uniform', 'normal', 'uniform', 'uniform']

p_sigma = [0., 0.5, 0., 0., 0., 0.]

Thep_formsets whether the prior probabilityp[i]is logarithmic (value = 10**p[i]) or linear (value = p[i]). p_formpriorsets the functional form of the prior proba-bility distribution. It should be set to normal oruniform for a Gaussian or rectan-gular distribution, respectively. p_sigma is the width of the Gaussian function.

• **SiDE** **on computer clusters.** SiDE has the option use_mpi where it is possible
to use MPI pools instead of Python threads. This alternative is useful for running
MCMC procedures on computer clusters using multiple nodes. To implement this,
ﬁrst it is necessary to add the runtime parameters as number of walkers (nwalkers)
and the number of walker steps in the main run (nsteps) in the ﬁtting Python script.

*# Runtime parameters*
nwalkers = 240
nburnin = 0