diff --git a/tardis/io/config_reader.py b/tardis/io/config_reader.py index 3be41bec68b..8ec20b5a6dd 100644 --- a/tardis/io/config_reader.py +++ b/tardis/io/config_reader.py @@ -30,9 +30,7 @@ def parse_convergence_section(convergence_section_dict): dictionary """ - - convergence_parameters = ['damping_constant', 'threshold', 'fraction', - 'hold_iterations'] + convergence_parameters = ['damping_constant', 'threshold'] for convergence_variable in ['t_inner', 't_rad', 'w']: if convergence_variable not in convergence_section_dict: @@ -260,9 +258,17 @@ def from_config_dict(cls, config_dict, validate=True, config_dirname=''): validated_config_dict['config_dirname'] = config_dirname montecarlo_section = validated_config_dict['montecarlo'] - montecarlo_section['convergence_strategy'] = ( - parse_convergence_section( - montecarlo_section['convergence_strategy'])) + if montecarlo_section['convergence_strategy']['type'] == "damped": + montecarlo_section['convergence_strategy'] = ( + parse_convergence_section( + montecarlo_section['convergence_strategy'])) + elif montecarlo_section['convergence_strategy']['type'] == "custom": + raise NotImplementedError( + 'convergence_strategy is set to "custom"; ' + 'you need to implement your specific convergence treatment') + else: + raise ValueError('convergence_strategy is not "damped" ' + 'or "custom"') return cls(validated_config_dict) diff --git a/tardis/io/schemas/montecarlo.yml b/tardis/io/schemas/montecarlo.yml index 7fed6027f9f..5aaf30cdbec 100644 --- a/tardis/io/schemas/montecarlo.yml +++ b/tardis/io/schemas/montecarlo.yml @@ -75,7 +75,7 @@ properties: convergence_strategy: oneOf: - $ref: '#/definitions/convergence_strategy/damped' - - $ref: '#/definitions/convergence_strategy/specific' + - $ref: '#/definitions/convergence_strategy/custom' default: type: 'damped' required: @@ -88,75 +88,17 @@ definitions: type: object additionalProperties: false properties: - type: - enum: - - damped - damping_constant: - type: number - default: 0.5 - description: damping constant - t_inner: - type: object - additionalProperties: false - properties: - damping_constant: - type: number - default: 0.5 - description: damping constant - threshold: - type: number - description: specifies the threshold that is taken as convergence (i.e. - 0.05 means that the value does not change more than 5%) - t_rad: - type: object - additionalProperties: false - properties: - damping_constant: - type: number - default: 0.5 - description: damping constant - threshold: - type: number - description: specifies the threshold that is taken as convergence (i.e. - 0.05 means that the value does not change more than 5%) - required: - - threshold - w: - type: object - additionalProperties: false - properties: - damping_constant: - type: number - default: 0.5 - description: damping constant - threshold: - type: number - description: specifies the threshold that is taken as convergence (i.e. - 0.05 means that the value does not change more than 5%) - required: - - threshold - lock_t_inner_cycles: - type: number - multipleOf: 1.0 - default: 1 - description: The number of cycles to lock the update of the inner boundary - temperature. This process helps with convergence. The default is to switch - it off (1 cycle) - t_inner_update_exponent: - type: number - default: -0.5 - description: L=4*pi*r**2*T^y - specific: - type: object - additionalProperties: false properties: type: enum: - - specific - threshold: - type: number - description: specifies the threshold that is taken as convergence (i.e. - 0.05 means that the value does not change more than 5%) + - damped + default: damped + stop_if_converged: + type: boolean + default: false + description: stop plasma iterations before number of specified + iterations are reached if the simulation is plasma and inner + boundary state is converged fraction: type: number default: 0.8 @@ -169,6 +111,15 @@ definitions: default: 3 description: the number of iterations that the convergence criteria need to be fulfilled before TARDIS accepts the simulation as converged + damping_constant: + type: number + default: 1.0 + description: damping constant + threshold: + type: number + default: 0.05 + description: specifies the threshold that is taken as convergence (i.e. + 0.05 means that the value does not change more than 5%) t_inner: type: object additionalProperties: false @@ -216,15 +167,21 @@ definitions: description: The number of cycles to lock the update of the inner boundary temperature. This process helps with convergence. The default is to switch it off (1 cycle) - damping_constant: - type: number - default: 0.5 - description: damping constant t_inner_update_exponent: type: number default: -0.5 description: L=4*pi*r**2*T^y - required: - - threshold - - fraction - - hold_iterations + custom: + type: object + additionalProperties: false + properties: + type: + enum: + - custom + description: Use this convergence_strategy for your specific needs. You + need to change the codebase accordingly + required: + - fraction + - damping_constant + - threshold + - hold_iterations diff --git a/tardis/io/tests/data/tardis_configv1_ascii_density_abund.yml b/tardis/io/tests/data/tardis_configv1_ascii_density_abund.yml index f922dbaf0a8..e0a6e9c684e 100644 --- a/tardis/io/tests/data/tardis_configv1_ascii_density_abund.yml +++ b/tardis/io/tests/data/tardis_configv1_ascii_density_abund.yml @@ -11,7 +11,7 @@ supernova: atom_data: kurucz_atom_pure_simple.h5 model: - + structure: type: file filename: density.dat @@ -35,7 +35,7 @@ montecarlo: seed: 23111963 no_of_packets : 1.0e+5 iterations: 20 - + black_body_sampling: start: 1 angstrom stop: 1000000 angstrom @@ -44,28 +44,15 @@ montecarlo: no_of_virtual_packets: 5 convergence_strategy: - type: specific + type: damped damping_constant: 1.0 threshold: 0.05 fraction: 0.8 hold_iterations: 3 t_inner: damping_constant: 1.0 - + spectrum: start : 500 angstrom stop : 20000 angstrom num: 10000 - - - - - - - - - - - - - diff --git a/tardis/io/tests/data/tardis_configv1_uniform_density.yml b/tardis/io/tests/data/tardis_configv1_uniform_density.yml index 675c3f2784b..a37b6ee8984 100644 --- a/tardis/io/tests/data/tardis_configv1_uniform_density.yml +++ b/tardis/io/tests/data/tardis_configv1_uniform_density.yml @@ -48,7 +48,7 @@ montecarlo: no_of_virtual_packets: 10 convergence_strategy: - type: specific + type: damped damping_constant: 1.0 threshold: 0.05 fraction: 0.8 diff --git a/tardis/io/tests/data/tardis_configv1_verysimple.yml b/tardis/io/tests/data/tardis_configv1_verysimple.yml index 7d1fea5e3ba..996f35434f6 100644 --- a/tardis/io/tests/data/tardis_configv1_verysimple.yml +++ b/tardis/io/tests/data/tardis_configv1_verysimple.yml @@ -39,6 +39,12 @@ montecarlo: iterations: 5 last_no_of_packets: 5.0e+5 no_of_virtual_packets: 5 + convergence_strategy: + type: damped + damping_constant: 0.5 + threshold: 0.05 + lock_t_inner_cycles: 1 + t_inner_update_exponent: -0.5 spectrum: start: 500 angstrom diff --git a/tardis/plasma/tests/data/plasma_test_config_lte.yml b/tardis/plasma/tests/data/plasma_test_config_lte.yml index 1e0d9e8b2af..c0bee588022 100644 --- a/tardis/plasma/tests/data/plasma_test_config_lte.yml +++ b/tardis/plasma/tests/data/plasma_test_config_lte.yml @@ -40,6 +40,12 @@ montecarlo: iterations: 10 last_no_of_packets: 1.0e+5 no_of_virtual_packets: 5 + convergence_strategy: + type: damped + damping_constant: 0.5 + threshold: 0.05 + lock_t_inner_cycles: 1 + t_inner_update_exponent: -0.5 spectrum: start: 500 angstrom diff --git a/tardis/plasma/tests/data/plasma_test_config_nlte.yml b/tardis/plasma/tests/data/plasma_test_config_nlte.yml index 35a844cd9cf..ab6417b08cb 100644 --- a/tardis/plasma/tests/data/plasma_test_config_nlte.yml +++ b/tardis/plasma/tests/data/plasma_test_config_nlte.yml @@ -42,6 +42,12 @@ montecarlo: iterations: 10 last_no_of_packets: 1.0e+5 no_of_virtual_packets: 5 + convergence_strategy: + type: damped + damping_constant: 0.5 + threshold: 0.05 + lock_t_inner_cycles: 1 + t_inner_update_exponent: -0.5 spectrum: start: 500 angstrom diff --git a/tardis/simulation/base.py b/tardis/simulation/base.py index 011e29ae9a9..8139af8c734 100644 --- a/tardis/simulation/base.py +++ b/tardis/simulation/base.py @@ -56,14 +56,19 @@ def __init__(self, iterations, model, plasma, runner, self.luminosity_nu_end = luminosity_nu_end self.luminosity_requested = luminosity_requested self.nthreads = nthreads - if convergence_strategy.type in ('damped', 'specific'): + if convergence_strategy.type in ('damped'): self.convergence_strategy = convergence_strategy self.converged = False self.consecutive_converges_count = 0 + elif convergence_strategy.type in ('custom'): + raise NotImplementedError( + 'Convergence strategy type is custom; ' + 'you need to implement your specific treatment!' + ) else: raise ValueError( 'Convergence strategy type is ' - 'neither damped nor specific ' + 'not damped or custom ' '- input is {0}'.format(convergence_strategy.type)) self._callbacks = OrderedDict() @@ -97,41 +102,37 @@ def _get_convergence_status(self, t_rad, w, t_inner, estimated_t_rad, convergence_t_inner = (abs(t_inner - estimated_t_inner) / estimated_t_inner).value - if self.convergence_strategy.type == 'specific': - fraction_t_rad_converged = ( - np.count_nonzero( - convergence_t_rad < self.convergence_strategy.t_rad.threshold) - / no_of_shells) - - t_rad_converged = ( - fraction_t_rad_converged > self.convergence_strategy.t_rad.threshold) - - fraction_w_converged = ( - np.count_nonzero( - convergence_w < self.convergence_strategy.w.threshold) - / no_of_shells) - - w_converged = ( - fraction_w_converged > self.convergence_strategy.w.threshold) - - t_inner_converged = ( - convergence_t_inner < self.convergence_strategy.t_inner.threshold) - - if np.all([t_rad_converged, w_converged, t_inner_converged]): - hold_iterations = self.convergence_strategy.hold_iterations - self.consecutive_converges_count += 1 - logger.info("Iteration converged {0:d}/{1:d} consecutive " - "times.".format(self.consecutive_converges_count, - hold_iterations + 1)) - # If an iteration has converged, require hold_iterations more - # iterations to converge before we conclude that the Simulation - # is converged. - return self.consecutive_converges_count == hold_iterations + 1 - else: - self.consecutive_converges_count = 0 - return False - + fraction_t_rad_converged = ( + np.count_nonzero( + convergence_t_rad < self.convergence_strategy.t_rad.threshold) + / no_of_shells) + + t_rad_converged = ( + fraction_t_rad_converged > self.convergence_strategy.fraction) + + fraction_w_converged = ( + np.count_nonzero( + convergence_w < self.convergence_strategy.w.threshold) + / no_of_shells) + + w_converged = ( + fraction_w_converged > self.convergence_strategy.fraction) + + t_inner_converged = ( + convergence_t_inner < self.convergence_strategy.t_inner.threshold) + + if np.all([t_rad_converged, w_converged, t_inner_converged]): + hold_iterations = self.convergence_strategy.hold_iterations + self.consecutive_converges_count += 1 + logger.info("Iteration converged {0:d}/{1:d} consecutive " + "times.".format(self.consecutive_converges_count, + hold_iterations + 1)) + # If an iteration has converged, require hold_iterations more + # iterations to converge before we conclude that the Simulation + # is converged. + return self.consecutive_converges_count == hold_iterations + 1 else: + self.consecutive_converges_count = 0 return False def advance_state(self): @@ -147,7 +148,8 @@ def advance_state(self): estimated_t_rad, estimated_w = ( self.runner.calculate_radiationfield_properties()) estimated_t_inner = self.estimate_t_inner( - self.model.t_inner, self.luminosity_requested) + self.model.t_inner, self.luminosity_requested, + t_inner_update_exponent=self.convergence_strategy.t_inner_update_exponent) converged = self._get_convergence_status(self.model.t_rad, self.model.w, @@ -163,9 +165,12 @@ def advance_state(self): self.convergence_strategy.t_rad.damping_constant) next_w = self.damped_converge( self.model.w, estimated_w, self.convergence_strategy.w.damping_constant) - next_t_inner = self.damped_converge( - self.model.t_inner, estimated_t_inner, - self.convergence_strategy.t_inner.damping_constant) + if (self.iterations_executed + 1) % self.convergence_strategy.lock_t_inner_cycles == 0: + next_t_inner = self.damped_converge( + self.model.t_inner, estimated_t_inner, + self.convergence_strategy.t_inner.damping_constant) + else: + next_t_inner = self.model.t_inner self.log_plasma_state(self.model.t_rad, self.model.w, self.model.t_inner, next_t_rad, next_w, @@ -211,10 +216,13 @@ def iterate(self, no_of_packets, no_of_virtual_packets=0, last_run=False): def run(self): start_time = time.time() - while self.iterations_executed < self.iterations-1 and not self.converged: + while self.iterations_executed < self.iterations-1: self.iterate(self.no_of_packets) self.converged = self.advance_state() self._call_back() + if self.converged: + if self.convergence_strategy.stop_if_converged: + break # Last iteration self.iterate(self.last_no_of_packets, self.no_of_virtual_packets, True)