Merge pull request #263 from Magritte-code/min_line_opacity

Changed default min_line_opacity parameter value

Author: ThomasCeulemans

docs/src/4_extra_options/NLTE_in_low_density.rst

@@ -3,14 +3,14 @@ NLTE in very low density regimes
 
 In low density regimes, the level populations might start to maser. 
 This means that the resulting line opacity (and optical depth) can become close to zero or even negative.
-Unfortunately, our solvers cannot handle this situation too well, and might produce unphysical results.
+Unfortunately, our solvers cannot handle a rapid change in source function (emissivity/opacity) too well, and might produce unphysical results.
 
 This is why we implement a workaround by setting the minimum allowed value for the line opacity (starting from version 0.7.0).
-This value is set to 1e-13 [W sr−1 m−3] by default, but can be changed by the user. Higher values might overestimate the impact of the less dense model regions on the intensity, while lower value can lead to numerical artifacts in the NLTE intensities.
+This value is set to 1e-10 [W sr−1 m−3] by default, but can be changed by the user. Higher values might overestimate the impact of the less dense model regions on the intensity, while lower value can lead to numerical artifacts in the NLTE intensities.
 
 .. code-block:: python
 
-    parameters.min_line_opacity = 1e-13#default
+    parameters.min_line_opacity = 1e-10#default
 
 In versions 0.5.3 until 0.6.0, the code implemented another workaround to avoid this situation by resetting specific levels at specific places to LTE.
 This might work in some cases (with the caveat of some artifacts being produced in the resulting spectrum), 

src/model/parameters/parameters.hpp

@@ -77,7 +77,7 @@ struct Parameters {
     Real min_rel_pop_for_convergence = 1.0e-10;
     Real pop_prec                    = 1.0e-6;
     Real min_opacity                 = 1.0e-26;
-    long double min_line_opacity     = 1.0e-13;
+    long double min_line_opacity     = 1.0e-10;
     Real min_dtau                    = 1.0e-15;
     Real population_inversion_fraction =
         -1.0; // previously 1.01; // threshold factor for population inversion required

tests/benchmarks/analytic/all_zero_single_ray.py

@@ -88,6 +88,8 @@ def run_model (nosave=False):
     timer1.start()
     model = magritte.Model (modelFile)
     timer1.stop()
+    
+    model.parameters.min_line_opacity = 1.0e-13
 
     timer2 = tools.Timer('setting model')
     timer2.start()
@@ -96,6 +98,7 @@ def run_model (nosave=False):
     model.compute_LTE_level_populations   ()
     timer2.stop()
 
+
     timer3 = tools.Timer('shortchar 0  ')
     timer3.start()
     model.compute_radiation_field_shortchar_order_0 ()

tests/benchmarks/analytic/density_distribution_1D.py

@@ -96,6 +96,8 @@ def run_model (a_or_b, nosave=False, use_widgets=True):
     model = magritte.Model (modelFile)
     timer1.stop()
 
+    model.parameters.min_line_opacity = 1.0e-13
+
     timer2 = tools.Timer('setting model')
     timer2.start()
     model.compute_spectral_discretisation ()

tests/benchmarks/analytic/density_distribution_1D_image.py

@@ -96,6 +96,8 @@ def run_model (a_or_b, nosave=False, use_widgets=True):
     model = magritte.Model (modelFile)
     timer1.stop()
 
+    model.parameters.min_line_opacity = 1.0e-13
+
     fcen = model.lines.lineProducingSpecies[0].linedata.frequency[0]
     dd = 3.0e+3 / magritte.CC
     fmin = fcen - fcen*dd

tests/benchmarks/analytic/density_distribution_single_ray.py

@@ -95,6 +95,8 @@ def run_model (a_or_b, nosave=False):
     model = magritte.Model (modelFile)
     timer1.stop()
 
+    model.parameters.min_line_opacity = 1.0e-13
+
     timer2 = tools.Timer('setting model')
     timer2.start()
     model.compute_spectral_discretisation ()

tests/benchmarks/analytic/overview.ipynb

@@ -13,7 +13,7 @@
 import magritte.setup    as setup
 import magritte.core     as magritte
 
-VERSION = "0.7.0"#reference version of Magritte to compare against; should correspond to a commited intensity file
+VERSION = "0.7.2"#reference version of Magritte to compare against; should correspond to a commited intensity file
 
 path = os.path.dirname(os.path.realpath(__file__))
 

tests/benchmarks/numeric/run_phantom_3D_model_reduced.py

@@ -103,6 +103,8 @@ def run_model (a_or_b, nosave=False):
     model.compute_LTE_level_populations   ()
     timer2.stop()
 
+    model.parameters.min_line_opacity = 1e-13
+
     timer3 = tools.Timer('running model')
     timer3.start()
     model.compute_level_populations (True, 100)

tests/benchmarks/numeric/vanZadelhoff_1_1D.py

@@ -103,6 +103,8 @@ def run_model (a_or_b, nosave=False):
     model.compute_LTE_level_populations   ()
     timer2.stop()
 
+    model.parameters.min_line_opacity = 1.0e-13
+
     timer3 = tools.Timer('running model')
     timer3.start()
     model.compute_level_populations_sparse (True, 100)